From 517b1918fb47e06d88b8f784e1d8f54f27561edb Mon Sep 17 00:00:00 2001
From: davidmcdonagh <davidnmcdonagh@gmail.com>
Date: Mon, 9 Sep 2024 14:06:25 +0100
Subject: [PATCH] Refactor asu reflection prediction for laue data

---
 .../integration/fit/tof_line_profile.py       |  3 +-
 .../integration/tof/tof_integration.h         | 32 ++++++-------------
 .../spot_prediction/reflection_predictor.h    | 10 +++---
 3 files changed, 16 insertions(+), 29 deletions(-)

diff --git a/src/dials/algorithms/integration/fit/tof_line_profile.py b/src/dials/algorithms/integration/fit/tof_line_profile.py
index b2613e5fb9..60410b053c 100644
--- a/src/dials/algorithms/integration/fit/tof_line_profile.py
+++ b/src/dials/algorithms/integration/fit/tof_line_profile.py
@@ -122,7 +122,8 @@ def compute_line_profile_data_for_reflection(
         l.fit()
         line_profile = l.result()
         fit_intensity = integrate.simpson(line_profile, x=tof)
-    except ValueError:
+    except ValueError as e:
+        print("fit error", e)
         return [], [], [], [], -1, -1, -1, -1
 
     if n_background > 0:
diff --git a/src/dials/algorithms/integration/tof/tof_integration.h b/src/dials/algorithms/integration/tof/tof_integration.h
index 71e23f654e..eda1bf7948 100644
--- a/src/dials/algorithms/integration/tof/tof_integration.h
+++ b/src/dials/algorithms/integration/tof/tof_integration.h
@@ -45,9 +45,9 @@ namespace dials { namespace algorithms {
   using scitbx::constants::Planck;
 
   void get_asu_reflections(af::shared<cctbx::miller::index<int> > indices,
-                           af::shared<cctbx::miller::index<int> > predicted_indices,
+                           af::shared<cctbx::miller::index<int> > asu_predicted_indices,
                            af::shared<double> wavelengths,
-                           af::shared<double> predicted_wavelengths,
+                           af::shared<double> asu_predicted_wavelengths,
                            af::shared<bool> asu_reflection,
                            cctbx::sgtbx::space_group space_group
 
@@ -58,7 +58,8 @@ namespace dials { namespace algorithms {
 
     DIALS_ASSERT(indices.size() == asu_reflection.size());
     DIALS_ASSERT(indices.size() == wavelengths.size());
-    DIALS_ASSERT(predicted_indices.size() == predicted_wavelengths.size());
+    DIALS_ASSERT(asu_predicted_indices.size() == asu_predicted_wavelengths.size());
+
     const char* hall_symbol = space_group.type().hall_symbol().c_str();
     const gemmi::SpaceGroup* gemmi_sg_ptr =
       gemmi::find_spacegroup_by_ops(gemmi::symops_from_hall(hall_symbol));
@@ -75,29 +76,14 @@ namespace dials { namespace algorithms {
       int isym = hkl_mover.move_to_asu(hkl);
       merged_hkls[i_refl] = cctbx::miller::index<>(hkl[0], hkl[1], hkl[3]);
     }
-    int matched = 0;
-    int unmatched = 0;
-    for (std::size_t i = 0; i < indices.size(); ++i) {
-      bool match = false;
-      for (std::size_t j = 0; j < predicted_indices.size(); ++j) {
-        if (indices[i] == predicted_indices[j]) {
-          match = true;
-          matched++;
-          break;
-        }
-      }
-      if (!match) {
-        unmatched++;
-      }
-    }
 
-    for (std::size_t i = 0; i < predicted_indices.size(); ++i) {
-      cctbx::miller::index<> p_hkl = predicted_indices[i];
+    for (std::size_t i = 0; i < asu_predicted_indices.size(); ++i) {
+      cctbx::miller::index<> p_hkl = asu_predicted_indices[i];
       int closest_match = -1;
       double min_wl_diff = -1;
-      for (std::size_t j = 0; j < indices.size(); ++j) {
-        if (p_hkl == indices[j]) {
-          double wl_diff = std::abs(wavelengths[j] - predicted_wavelengths[j]);
+      for (std::size_t j = 0; j < merged_hkls.size(); ++j) {
+        if (p_hkl == merged_hkls[j]) {
+          double wl_diff = std::abs(wavelengths[j] - asu_predicted_wavelengths[i]);
           if (min_wl_diff < 0 || wl_diff < min_wl_diff) {
             closest_match = j;
           }
diff --git a/src/dials/algorithms/spot_prediction/reflection_predictor.h b/src/dials/algorithms/spot_prediction/reflection_predictor.h
index 83d9e615b5..f34a167a7d 100644
--- a/src/dials/algorithms/spot_prediction/reflection_predictor.h
+++ b/src/dials/algorithms/spot_prediction/reflection_predictor.h
@@ -1372,8 +1372,6 @@ namespace dials { namespace algorithms {
       cctbx::miller::index_generator indices =
         cctbx::miller::index_generator(unit_cell_, space_group_type_, false, dmin_);
 
-      af::shared<miller_index> indices_arr = indices.to_array();
-
       af::reflection_table table;
       af::shared<double> wavelength_column;
       table["wavelength_cal"] = wavelength_column;
@@ -1381,9 +1379,11 @@ namespace dials { namespace algorithms {
       table["s0_cal"] = s0_column;
       laue_prediction_data predictions(table);
 
-      for (std::size_t i = 0; i < indices_arr.size(); ++i) {
-        miller_index h = indices_arr[i];
-
+      for (;;) {
+        miller_index h = indices.next();
+        if (h.is_zero()) {
+          break;
+        }
         vec3<double> q = setting_rotation * rotation * fixed_rotation * ub_ * h;
 
         // Calculate the wavelength required to meet the diffraction condition