rugnux: a rotation frame the sweep's cell refuses is tried on that cell scaled

A rotation frame is refined on the sweep's one cell (orientation and beam per frame, never the cell), and a crystal whose cell grows with the dose leaves that cell behind: the per-frame gate (20 % of the in-resolution, non-ice spots on the lattice) then refuses the frames it misses as if the crystal had gone. A frame that fails the gate is now tried once more on the sweep's cell scaled isotropically - steps of 0.1 %, out to +-2 %, nearest first - to the scale that puts most of its non-ice spots on the lattice, re-refined there, and integrated on it where it then clears the gate. The scale must gain more spots than sqrt(count), so forty trial scales cannot turn a chance spot on a sparse frame into a frame. Ice-ring spots are not counted, as the gate does not count them. A frame the sweep's cell already admits is untouched, and only the per-image integration path does this (not the first-pass validation scoring). Found on 5ky6 (P 1 21 1, 564 x 0.5 deg, B rising by 46 A^2 over the sweep), which part B (5b9144410, decisions before pass 2 from the first part of the sweep) turned from P 1 21 1 into P 1 2 1. The mechanism: - The first pass now indexes on [0, 141 deg). A per-frame probe of the isotropic scale that puts most spots on the canonical lattice reads -0.1 % on frame 0, +0.4 % at 50 deg, >= +0.6 % from 115 deg on (at frame 230: 62 spots on the lattice at the sweep's cell, 112 at +0.6 %). Sampling the whole sweep (before part B) handed pass 2 a cell averaged over the dose, so frames to ~147 deg and part of the last wedge still cleared the gate; the early-range cell is the low-dose cell, and every frame from 113 deg on was refused (221 of 564 frames merged, completeness 83 %). - The axial 0k0 row is recorded in those frames: the merge held 1 screw-absent axial reflection instead of 15, so the 2_1 screw could not be decided (P2 by default), and R-free 0.243 -> 0.263. - Separately, on the early range the FFT candidate matched a false C-orthorhombic setting (gamma 92 deg) whose refutation fell straight to triclinic instead of P-monoclinic, so the cell was refined free (cell dev 0.64 -> 1.03 %). Not changed here (see the report). Integrating every frame from the sweep's cell instead (the sparse-lattice switch forced on) also recovers P 1 21 1, but at 1.63 A with ISa 5.0 and R_meas 45 %: the late frames are then integrated on a cell they have left. Result (GPU, battery verdicts against the rc175 full run 7c9a6c7d1): 5ky6 fail P 1 2 1, d 1.50, R-free 0.263, compl 83.0 % -> pass P 1 21 1, d 1.62, R-free 0.234, compl 99.6 % (pre-part-B: P 1 21 1, 1.54, 0.244, 96.9 %) 9rci pass -> pass, compl 86.5 -> 94.7 %, d 1.81 -> 1.76 9qw8 pass -> pass, R-free 0.329 -> 0.326, ISa 5.2 -> 5.5 8v4o fail -> fail (its part-B lattice loss is untouched), compl 22.3 -> 23.1 % p.mtz md5 unchanged: 5nw5, insu_I_x06da_low_isa, nothing_1, nothing_2, myob_x06da_powder_2, 7atg, 5reo, 6fvz, 6hwj, 6pxc, 8tha, 9crw, 9jzo, 6zqr, 6rym, thau_x10sa_0p1deg, cytc_x10sa. Catch2: Rugnux_Rotation, RotationIndexer* pass. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01SVmAWnzCmRKAXVUCdc4iNi
2026-10-08 21:56:46 +02:00
parent 7a3bf6c5cb
commit a825fb6cba
3 changed files with 66 additions and 8 deletions
@@ -39,7 +39,7 @@ A spot is indexed if $\delta^2 < \tau^2$, where $\tau$ is the configured toleran

 For indexed spots, the reciprocal lattice point $\mathbf{p} = h\mathbf{a}^*+k\mathbf{b}^*+l\mathbf{c}^*$ is used to compute $\Delta_\mathrm{Ewald}(\mathbf{p})$ (stored as a diagnostic and later used in profile-radius estimation).

-A frame is taken to be this crystal's when at least a fraction $g = 0.20$ of its in-resolution, non-ice spots index. On rotation data that decision is what admits the frame to integration, so its denominator matters: every spot handed to it that is not a reflection of this crystal argues against the frame. Where it cannot do that job — a lattice whose pooled spots fall below $g$ and which fewer than half the validation frames clear — every frame is integrated instead: the frames that clear $g$ are then only the upper tail of the same sparse population, not the frames the crystal was in, and admitting them alone dropped most of a small-molecule sweep and its completeness with it.
+A frame is taken to be this crystal's when at least a fraction $g = 0.20$ of its in-resolution, non-ice spots index. On rotation data that decision is what admits the frame to integration, so its denominator matters: every spot handed to it that is not a reflection of this crystal argues against the frame. Where it cannot do that job — a lattice whose pooled spots fall below $g$ and which fewer than half the validation frames clear — every frame is integrated instead: the frames that clear $g$ are then only the upper tail of the same sparse population, not the frames the crystal was in, and admitting them alone dropped most of a small-molecule sweep and its completeness with it. A rotation frame is refined on the sweep's one cell, and a crystal whose cell grows with the dose leaves that cell behind, so a frame that fails $g$ is tried once more on the sweep's cell scaled isotropically (steps of 0.1 %, out to ±2 %) to the scale that puts most of its non-ice spots on the lattice, and is integrated on that cell where it then clears $g$. The scale has to gain more spots than the count's own noise, $\sqrt{n}$, so forty trial scales cannot turn a chance spot into a frame. A frame the sweep's cell already admits is integrated exactly as before. Measured on a $P2_1$ crystal losing 46 Å² over the sweep, the best scale walks from −0.1 % on the first frame to beyond +0.6 % at 115°, and without it the frames past that point — and the axial row that decides its screw axis — were refused.

 That test decides a **frame**. Whether a rotation run has a lattice **at all** is decided on the spots instead. A frame count comes from serial crystallography, where each image is its own experiment; a rotation sweep is one crystal and one orientation matrix, its frames are not independent of each other, and what such a count mostly measures is how many spots happen to land on a frame — a sweep carrying four spots an image cannot reach a six-spot bar on three frames in four however right the lattice is. The refusal therefore compares the fraction of *all* spots in the sampled frames that the lattice explains against what the same lattice explains when each frame's spots are put at **another frame's angle**: same lattice, same spots, same detector, same refinement, with only the claim that these spots were seen at *these* angles removed. That difference is the evidence, and it carries no spots-per-frame number anywhere, so nothing has to be chosen for a crystal that diffracts weakly. Measured over a hundred datasets the permuted level never exceeds 2.6 % and the smallest true margin is seventeen points. `--min-indexed-spots` (default 6, floor 4 — four is where a lattice stops being fitted by any three spots) still sets the reported indexing rate and the count the first pass *ranks* candidate lattices by; every rescue and every arbiter still counts frames.

@@ -698,7 +698,7 @@ void IndexAndRefine::ProbeSupercellFrame(const DataMessage &msg, BraggPrediction

 std::optional<IndexAndRefine::IndexingOutcome>
 IndexAndRefine::DetermineRefineAnalyze(DataMessage &msg, const SpotFindingSettings &spot_finding_settings,
-                                      int64_t *spots_on_lattice) {
+                                      int64_t *spots_on_lattice, bool follow_cell_drift) {
    if (!indexer_ || !spot_finding_settings.indexing)
        return std::nullopt;

@@ -712,7 +712,12 @@ IndexAndRefine::DetermineRefineAnalyze(DataMessage &msg, const SpotFindingSettin
    if (!outcome.lattice_candidate)
        return std::nullopt;

-    if (experiment.GetIndexingSettings().GetGeomRefinementAlgorithm() != GeomRefinementAlgorithmEnum::None)
+    // Kept for the cell-drift retry below, on the only path that makes it.
+    std::optional<IndexingOutcome> unrefined;
+    if (rotation_indexer && follow_cell_drift)
+        unrefined = outcome;
+    const bool refine = experiment.GetIndexingSettings().GetGeomRefinementAlgorithm() != GeomRefinementAlgorithmEnum::None;
+    if (refine)
        RefineGeometryIfNeeded(msg, outcome);

    if (!outcome.lattice_candidate.has_value())
@@ -721,8 +726,59 @@ IndexAndRefine::DetermineRefineAnalyze(DataMessage &msg, const SpotFindingSettin
    // AnalyzeIndexing answers "is this frame worth integrating"; msg.indexing_result carries the
    // stricter "does this frame index on its own", which on rotation is not the same question.
    if (!AnalyzeIndexing(msg, outcome.experiment, *outcome.lattice_candidate, outcome.extra_lattice_candidates,
-                         spots_on_lattice, integrate_every_frame_))
-        return std::nullopt;
+                         spots_on_lattice, integrate_every_frame_)) {
+        // A rotation frame is refined on the sweep's cell, and a crystal whose cell grows with the dose
+        // leaves that cell behind: the frames it misses are refused as if the crystal had gone. So the
+        // other hypothesis - this frame's cell is the sweep's, scaled - is tried on a frame that is
+        // refused, and the frame is integrated on it where it then indexes. Measured on a P2_1 crystal
+        // losing B by 46 A^2 over the sweep: the scale that puts most spots on the lattice walks from
+        // -0.1 % on the first frame to +0.4 % at 50 deg and past +0.6 % at 115 deg, where the sweep's
+        // cell puts 60-80 of ~400 spots on the lattice and that one 110-150, so the frames past it were
+        // refused - and with them the axial row that decides the screw axis. A frame the sweep's cell
+        // already integrates is left exactly as it was.
+        if (!unrefined)
+            return std::nullopt;
+        const float tol = experiment.GetIndexingSettings().GetTolerance();
+        // The spots the gate counts: an ice ring sits at one resolution, and a cell scaled to put its
+        // nodes on the ring would otherwise be scored on it.
+        const bool index_ice_rings = experiment.GetIndexingSettings().GetIndexIceRings();
+        std::vector<SpotToSave> counted;
+        for (const auto &sp : msg.spots)
+            if (index_ice_rings || !sp.ice_ring)
+                counted.push_back(sp);
+        const CrystalLattice &refined = *outcome.lattice_candidate;
+        const auto geom = outcome.experiment.GetDiffractionGeometry();
+        auto scaled = [](const CrystalLattice &l, float s) {
+            return CrystalLattice(l.Vec0() * s, l.Vec1() * s, l.Vec2() * s);
+        };
+        // Steps of 0.1 %, out to +-2 %, nearest first, so a tie goes to the smaller change. The scale
+        // has to put more spots on the lattice than the count's own noise, sqrt(count), or forty trial
+        // scales would turn one chance spot on a sparse frame into a frame that passes.
+        const int sweep_count = CountIndexedSpots(geom, refined, counted, tol * tol);
+        float best_scale = 1.0f;
+        int best_count = sweep_count;
+        for (int i = 1; i <= 20; i++)
+            for (const int sign : {1, -1}) {
+                const float s = 1.0f + 0.001f * static_cast<float>(sign * i);
+                const int n = CountIndexedSpots(geom, scaled(refined, s), counted, tol * tol);
+                if (n > best_count) {
+                    best_count = n;
+                    best_scale = s;
+                }
+            }
+        if (best_count <= sweep_count + std::sqrt(static_cast<float>(best_count)))
+            return std::nullopt;
+        outcome = *unrefined;
+        outcome.lattice_candidate = scaled(*unrefined->lattice_candidate, best_scale);
+        for (auto &el : outcome.extra_lattice_candidates)
+            el = scaled(el, best_scale);
+        if (refine)
+            RefineGeometryIfNeeded(msg, outcome);
+        if (!outcome.lattice_candidate.has_value()
+            || !AnalyzeIndexing(msg, outcome.experiment, *outcome.lattice_candidate,
+                                outcome.extra_lattice_candidates, spots_on_lattice, integrate_every_frame_))
+            return std::nullopt;
+    }

    {
        std::unique_lock ul(reflections_mutex);
@@ -737,7 +793,7 @@ void IndexAndRefine::ProcessImage(DataMessage &msg,
                                  BraggPrediction &prediction,
                                  const BraggIntegrateFn &integrate,
                                  const BraggIntegrateFn &probe_integrate) {
-    auto outcome = DetermineRefineAnalyze(msg, spot_finding_settings);
+    auto outcome = DetermineRefineAnalyze(msg, spot_finding_settings, nullptr, /*follow_cell_drift=*/true);
    if (outcome && spot_finding_settings.quick_integration)
        QuickPredictAndIntegrate(msg, spot_finding_settings, prediction, integrate, *outcome, probe_integrate);
 }
@@ -139,10 +139,12 @@ class IndexAndRefine {
    // Shared indexing path: determine the lattice/symmetry, refine geometry, and run AnalyzeIndexing.
    // Returns the outcome (ready for integration) when the frame indexes, nullopt otherwise. Both the
    // real per-image ProcessImage and the first-pass scheme validation go through this, so they cannot
-    // diverge.
+    // diverge. follow_cell_drift (ProcessImage only): a rotation frame the sweep's cell refuses is tried
+    // again on that cell scaled to the frame's own spots - the crystal's cell grows with the dose.
    std::optional<IndexingOutcome> DetermineRefineAnalyze(DataMessage &msg,
                                                          const SpotFindingSettings &spot_finding_settings,
-                                                          int64_t *spots_on_lattice = nullptr);
+                                                          int64_t *spots_on_lattice = nullptr,
+                                                          bool follow_cell_drift = false);
    void RefineGeometryIfNeeded(DataMessage &msg, IndexingOutcome &outcome);
    void QuickPredictAndIntegrate(DataMessage &msg,
                                 const SpotFindingSettings &spot_finding_settings,