indexamajig: Add --wavelength-estimate

This is a more sensible non-indexer-specific and non-Xray-specific
replacement for --pinkIndexer-override-photon-energy

4 files changed, 39 insertions, 28 deletions

diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c
index 094ed3d7..c06dda26 100644
--- a/libcrystfel/src/index.c
+++ b/libcrystfel/src/index.c
@@ -68,6 +68,7 @@ struct _indexingprivate
 	IndexingFlags flags;
 	UnitCell *target_cell;
 	double tolerance[6];
+	double wavelength_estimate;
 
 	int n_methods;
 	IndexingMethod *methods;
@@ -234,6 +235,7 @@ static char *friendly_indexer_name(IndexingMethod m)
 
 static void *prepare_method(IndexingMethod *m, UnitCell *cell,
                             const DataTemplate *dtempl,
+                            double wavelength_estimate,
                             struct xgandalf_options *xgandalf_opts,
                             struct pinkIndexer_options* pinkIndexer_opts,
                             struct felix_options *felix_opts,
@@ -283,7 +285,7 @@ static void *prepare_method(IndexingMethod *m, UnitCell *cell,
 
 		case INDEXING_PINKINDEXER :
 		priv = pinkIndexer_prepare(m, cell, pinkIndexer_opts,
-		                           dtempl);
+		                           dtempl, wavelength_estimate);
 		break;
 
 		default :
@@ -358,6 +360,7 @@ IndexingMethod *parse_indexing_methods(const char *method_list,
 IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell,
                                 const DataTemplate *dtempl,
                                 float *tols, IndexingFlags flags,
+                                double wavelength_estimate,
                                 struct taketwo_options *ttopts,
                                 struct xgandalf_options *xgandalf_opts,
                                 struct pinkIndexer_options *pinkIndexer_opts,
@@ -427,6 +430,7 @@ IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell,
 
 		ipriv->engine_private[i] = prepare_method(&methods[i], cell,
 		                                          dtempl,
+		                                          wavelength_estimate,
 		                                          xgandalf_opts,
 		                                          pinkIndexer_opts,
 		                                          felix_opts,
@@ -469,6 +473,7 @@ IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell,
 	ipriv->methods = methods;
 	ipriv->n_methods = n;
 	ipriv->flags = flags;
+	ipriv->wavelength_estimate = wavelength_estimate;
 
 	if ( cell != NULL ) {
 		ipriv->target_cell = cell_new_from_cell(cell);
@@ -940,6 +945,18 @@ void index_pattern_3(struct image *image, IndexingPrivate *ipriv, int *ping,
 	image->crystals = NULL;
 	image->n_crystals = 0;
 
+	if ( !isnan(ipriv->wavelength_estimate) ) {
+		if ( !within_tolerance(image->lambda,
+		                          ipriv->wavelength_estimate,
+		                          10.0) )
+		{
+			   ERROR("WARNING: Wavelength for %s %s (%e) differs by "
+			         "more than 10%% from estimated value (%e)\n",
+			         image->filename, image->ev,
+			         image->lambda, ipriv->wavelength_estimate);
+		}
+	}
+
 	orig = image->features;
 
 	for ( n=0; n<ipriv->n_methods; n++ ) {
diff --git a/libcrystfel/src/index.h b/libcrystfel/src/index.h
index 3edfc4cd..708c04a4 100644
--- a/libcrystfel/src/index.h
+++ b/libcrystfel/src/index.h
@@ -9,7 +9,7 @@
  * Copyright © 2012 Lorenzo Galli
  *
  * Authors:
- *   2010-2020 Thomas White <taw@physics.org>
+ *   2010-2021 Thomas White <taw@physics.org>
  *   2010      Richard Kirian
  *   2012      Lorenzo Galli
  *   2015      Kenneth Beyerlein <kenneth.beyerlein@desy.de>
@@ -168,6 +168,7 @@ extern IndexingPrivate *setup_indexing(const char *methods,
                                        const DataTemplate *dtempl,
                                        float *ltl,
                                        IndexingFlags flags,
+                                       double wavelength_estimate,
                                        struct taketwo_options *ttopts,
                                        struct xgandalf_options *xgandalf_opts,
                                        struct pinkIndexer_options *pinkIndexer_opts,
diff --git a/libcrystfel/src/indexers/pinkindexer.c b/libcrystfel/src/indexers/pinkindexer.c
index ac4a42d2..ae4e5642 100644
--- a/libcrystfel/src/indexers/pinkindexer.c
+++ b/libcrystfel/src/indexers/pinkindexer.c
@@ -3,11 +3,12 @@
  *
  * Interface to PinkIndexer
  *
- * Copyright © 2017-2020 Deutsches Elektronen-Synchrotron DESY,
+ * Copyright © 2017-2021 Deutsches Elektronen-Synchrotron DESY,
  *                       a research centre of the Helmholtz Association.
  *
  * Authors:
  *   2017-2019 Yaroslav Gevorkov <yaroslav.gevorkov@desy.de>
+ *   2021 Thomas White <thomas.white@desy.de>
  *
  * This file is part of CrystFEL.
  *
@@ -52,7 +53,6 @@ struct pinkIndexer_options {
 	int min_peaks;
 	int no_check_indexed;
 	float reflectionRadius; /* In m^-1 */
-	float customPhotonEnergy;
 	float customBandwidth;
 	float maxRefinementDisbalance;
 };
@@ -96,6 +96,7 @@ int run_pinkIndexer(struct image *image, void *ipriv)
 	reciprocalPeaks_1_per_A_t* reciprocalPeaks_1_per_A = &(pinkIndexer_private_data->reciprocalPeaks_1_per_A);
 	float *intensities = pinkIndexer_private_data->intensities;
 
+	/* FIXME: Check if wavelength is too far from estimate */
 	int peakCountMax = image_feature_count(image->features);
 	if (peakCountMax < 5) {
 		int goodLatticesCount = 0;
@@ -195,17 +196,17 @@ static int want_center_adjustment(struct pinkIndexer_options *pinkIndexer_opts)
 
 void *pinkIndexer_prepare(IndexingMethod *indm, UnitCell *cell,
                           struct pinkIndexer_options *pinkIndexer_opts,
-                          const DataTemplate *dtempl)
+                          const DataTemplate *dtempl,
+                          double wavelength_estimate)
 {
-	if ( !data_template_has_fixed_wavelength(dtempl)
-	     && (pinkIndexer_opts->customPhotonEnergy <= 0.0) )
-	{
-		ERROR("Geometry file refers to image metadata for "
-		      "wavelength.\n");
-		ERROR("To use PinkIndexer, specify a constant "
-		      "wavelength in the geometry file, or use "
-		      "--pinkIndexer-override-photon-energy.\n");
+	float beamEenergy_eV;
+
+	if ( isnan(wavelength_estimate) ) {
+		ERROR("PinkIndexer requires a wavelength estimate.  "
+		      "Try again with --wavelength-estimate=xx\n");
 		return NULL;
+	} else {
+		beamEenergy_eV = J_to_eV(ph_lambda_to_en(wavelength_estimate));
 	}
 
 	if ( !data_template_has_fixed_geometry(dtempl)
@@ -261,11 +262,7 @@ void *pinkIndexer_prepare(IndexingMethod *indm, UnitCell *cell,
 	}
 
 	/* FIXME: Beam gone */
-	float beamEenergy_eV = beam->photon_energy;
 	float nonMonochromaticity = beam->bandwidth*5;
-	if(pinkIndexer_opts->customPhotonEnergy > 0){
-		beamEenergy_eV = pinkIndexer_opts->customPhotonEnergy;
-	}
 	if(pinkIndexer_opts->customBandwidth >= 0){
 		nonMonochromaticity = pinkIndexer_opts->customBandwidth;
 	}
@@ -462,8 +459,6 @@ static void pinkIndexer_show_help()
 "     --pinkIndexer-max-refinement-disbalance=n\n"
 "                           Maximum disbalance after refinement:\n"
 "                            0 (no disbalance) to 2 (extreme disbalance), default 0.4\n"
-"     --pinkIndexer-override-photon-energy=ev\n"
-"                           Mean energy in eV to use for indexing.\n"
 "     --pinkIndexer-override-bandwidth=n\n"
 "                           Bandwidth in (delta energy)/(mean energy) to use for indexing.\n"
 "     --pinkIndexer-override-visible-energy-range=min-max\n"
@@ -492,7 +487,6 @@ int pinkIndexer_default_options(PinkIndexerOptions **opts_ptr)
 	opts->no_check_indexed = 0;
 	opts->min_peaks = 2;
 	opts->reflectionRadius = -1;
-	opts->customPhotonEnergy = -1;
 	opts->customBandwidth = -1;
 	opts->maxRefinementDisbalance = 0.4;
 
@@ -588,12 +582,9 @@ static error_t pinkindexer_parse_arg(int key, char *arg,
 		break;
 
 		case 11 :
-		if (sscanf(arg, "%f", &(*opts_ptr)->customPhotonEnergy) != 1)
-		{
-			ERROR("Invalid value for --pinkIndexer-override-photon-energy\n");
-			return EINVAL;
-		}
-		break;
+		ERROR("Please use --wavelength-estimate instead of "
+		      "--pinkIndexer-override-photon-energy.\n");
+		return EINVAL;
 
 		case 12 :
 		if (sscanf(arg, "%f", &(*opts_ptr)->customBandwidth) != 1)
diff --git a/libcrystfel/src/indexers/pinkindexer.h b/libcrystfel/src/indexers/pinkindexer.h
index 57df5f62..fc2714c8 100644
--- a/libcrystfel/src/indexers/pinkindexer.h
+++ b/libcrystfel/src/indexers/pinkindexer.h
@@ -3,11 +3,12 @@
  *
  * Interface to PinkIndexer
  *
- * Copyright © 2017-2020 Deutsches Elektronen-Synchrotron DESY,
+ * Copyright © 2017-2021 Deutsches Elektronen-Synchrotron DESY,
  *                       a research centre of the Helmholtz Association.
  *
  * Authors:
  *   2017-2019 Yaroslav Gevorkov <yaroslav.gevorkov@desy.de>
+ *   2021 Thomas White <thomas.white@desy.de>
  *
  * This file is part of CrystFEL.
  *
@@ -40,7 +41,8 @@ extern int run_pinkIndexer(struct image *image, void *ipriv);
 
 extern void *pinkIndexer_prepare(IndexingMethod *indm, UnitCell *cell,
                                  struct pinkIndexer_options *pinkIndexer_opts,
-                                 const DataTemplate *dtempl);
+                                 const DataTemplate *dtempl,
+                                 double wavelength_estimate);
 
 extern void pinkIndexer_cleanup(void *pp);