From e6bbfea4c4daaba9aa2183499e8a186eec2cd02c Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Wed, 7 Oct 2026 13:40:41 +0200 Subject: [PATCH] Indexing: a coplanar cell leaves the frame unindexed instead of throwing The online broker cancelled a collection with "Crystal lattice is coplanar and has no reciprocal cell". The candidate filter in PostIndexingRefinement (used by FFT, FFTW and ffbidx) tests lengths, angles and spot count only; three rows 120 deg apart in one plane pass all of them, and spots along the plane normal index as (0,0,0). Such a candidate - from ffbidx, or one the least-squares refinement walked flat - went on to LatticeSearch and XtalOptimizer, which refuses it and leaves it in place, and then reached Astar() in AnalyzeIndexing/prediction, which throws; the receiver turns that into a cancelled acquisition. - Refine() rejects a candidate below MIN_BASIS_VOLUME_FRACTION, like any other bad candidate. - XtalOptimizerInternal builds the refined lattice first and counts a coplanar result as a failed refinement, writing nothing back (the residual's volume clamp lets a solve end there and report success). Test: PostIndexingRefinement_CoplanarCandidateIsRejected threw the exact message before the fix. rugnux p.mtz unchanged on the three in-house reference sets at --prepass-fraction 1. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01SVmAWnzCmRKAXVUCdc4iNi --- .../geom_refinement/XtalOptimizer.cpp | 44 +++++++++++-------- .../indexing/PostIndexingRefinement.cpp | 6 +++ tests/IndexingUnitTest.cpp | 40 +++++++++++++++++ 3 files changed, 72 insertions(+), 18 deletions(-) diff --git a/image_analysis/geom_refinement/XtalOptimizer.cpp b/image_analysis/geom_refinement/XtalOptimizer.cpp index 23e671a2f..625430810 100644 --- a/image_analysis/geom_refinement/XtalOptimizer.cpp +++ b/image_analysis/geom_refinement/XtalOptimizer.cpp @@ -457,6 +457,31 @@ bool XtalOptimizerInternal(XtalOptimizerData &data, if (!summary.IsSolutionUsable()) return false; + // A solve can walk the basis flat and still report success - the volume clamp in XtalResidual + // keeps every step finite on the way. A coplanar cell has no reciprocal basis, so that is a + // failed refinement: refuse it like any other, before anything is written back. + CrystalLattice latt; + if (data.crystal_system == gemmi::CrystalSystem::Orthorhombic) + latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, PI / 2.0, PI / 2.0, PI / 2.0); + else if (data.crystal_system == gemmi::CrystalSystem::Tetragonal) { + latt_vec1[1] = latt_vec1[0]; + latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, PI / 2.0, PI / 2.0, PI / 2.0); + } else if (data.crystal_system == gemmi::CrystalSystem::Cubic) { + latt_vec1[1] = latt_vec1[0]; + latt_vec1[2] = latt_vec1[0]; + latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, PI / 2.0, PI / 2.0, PI / 2.0); + } else if (data.crystal_system == gemmi::CrystalSystem::Hexagonal) { + latt_vec1[1] = latt_vec1[0]; + latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1,PI / 2.0, PI / 2.0, 2.0 * PI / 3.0); + } else if (data.crystal_system == gemmi::CrystalSystem::Monoclinic) { + latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, PI / 2.0, latt_vec2[0], PI / 2.0); + } else { + // Triclinic via the same generic builder + latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, latt_vec2[0], latt_vec2[1], latt_vec2[2]); + } + if (latt.VolumeFraction() < MIN_BASIS_VOLUME_FRACTION) + return false; + if (data.refine_beam_center) { data.beam_corr_x = data.geom.GetBeamX_pxl() - beam[0]; data.beam_corr_y = data.geom.GetBeamY_pxl() - beam[1]; @@ -470,24 +495,7 @@ bool XtalOptimizerInternal(XtalOptimizerData &data, if (data.axis && data.refine_rotation_axis) data.axis.value().Axis(Coord(rot_vec[0], rot_vec[1], rot_vec[2])); - if (data.crystal_system == gemmi::CrystalSystem::Orthorhombic) - data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, PI / 2.0, PI / 2.0, PI / 2.0); - else if (data.crystal_system == gemmi::CrystalSystem::Tetragonal) { - latt_vec1[1] = latt_vec1[0]; - data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, PI / 2.0, PI / 2.0, PI / 2.0); - } else if (data.crystal_system == gemmi::CrystalSystem::Cubic) { - latt_vec1[1] = latt_vec1[0]; - latt_vec1[2] = latt_vec1[0]; - data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, PI / 2.0, PI / 2.0, PI / 2.0); - } else if (data.crystal_system == gemmi::CrystalSystem::Hexagonal) { - latt_vec1[1] = latt_vec1[0]; - data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1,PI / 2.0, PI / 2.0, 2.0 * PI / 3.0); - } else if (data.crystal_system == gemmi::CrystalSystem::Monoclinic) { - data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, PI / 2.0, latt_vec2[0], PI / 2.0); - } else { - // Triclinic via the same generic builder - data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, latt_vec2[0], latt_vec2[1], latt_vec2[2]); - } + data.latt = latt; return true; } catch (...) { // Convergence problems, likely not updated diff --git a/image_analysis/indexing/PostIndexingRefinement.cpp b/image_analysis/indexing/PostIndexingRefinement.cpp index a31bb7de3..ea906ae12 100644 --- a/image_analysis/indexing/PostIndexingRefinement.cpp +++ b/image_analysis/indexing/PostIndexingRefinement.cpp @@ -270,6 +270,12 @@ std::vector Refine(const std::vector &in_spots, cos_gamma > cos_min_angle || cos_gamma < cos_max_angle) continue; + // Three coplanar rows pass both tests above - 120 deg apart in one plane is inside any angle + // window - and the refinement can walk a candidate there. Such a cell has no reciprocal basis, + // so it is not a lattice: the same volume test as CrystalLattice::VolumeFraction(). + if (std::abs(cell_rows.determinant()) < MIN_BASIS_VOLUME_FRACTION * row_norms.prod()) + continue; + int64_t indexed_spot_count = 0; auto indexed_mask = ComputeIndexedMask(spots.topRows(nspots), cell_cols, p.indexing_tolerance, indexed_spot_count); diff --git a/tests/IndexingUnitTest.cpp b/tests/IndexingUnitTest.cpp index 33f43a920..2c05a0c01 100644 --- a/tests/IndexingUnitTest.cpp +++ b/tests/IndexingUnitTest.cpp @@ -785,3 +785,43 @@ TEST_CASE("FitSupercellProbe", "[SupercellProbe]") { // Too few reflections to fit a line: nothing. CHECK(FitSupercellProbe(SupercellProbeClass{}).b == 0.0); } + +TEST_CASE("PostIndexingRefinement_CoplanarCandidateIsRejected", "[Indexing]") { + // A candidate whose three rows lie in one plane passes every length and angle test - three rows + // 120 deg apart in a plane meet the 30-150 deg window - and every spot along the plane normal + // indexes as (0,0,0), so the spot count passes too. Such a cell has no reciprocal basis; handed + // on, it made the online analysis throw "Crystal lattice is coplanar" for the frame instead of + // leaving the frame unindexed. + const Coord a(50, 0, 0); + const Coord b(-25, 25.0f * std::sqrt(3.0f), 0); + const Coord c(-25, -25.0f * std::sqrt(3.0f), 0); + REQUIRE(CrystalLattice(a, b, c).VolumeFraction() < MIN_BASIS_VOLUME_FRACTION); + + std::vector spots; + for (int i = 1; i <= 40; i++) + spots.emplace_back(0.0f, 0.0f, 0.01f * static_cast(i)); + + Eigen::MatrixX3 cells(3, 3); + cells << a.x, a.y, a.z, + b.x, b.y, b.z, + c.x, c.y, c.z; + // A score below the refinement's distance floor, so the candidate reaches the filters as given. + Eigen::VectorX scores(1); + scores << -10.99f; + + RefineParameters parameters{ + .viable_cell_min_spots = 9, + .dist_tolerance_vs_reference = 0.05, + .reference_unit_cell = std::nullopt, + .min_length_A = 5, + .max_length_A = 500, + .min_angle_deg = 30, + .max_angle_deg = 150, + .indexing_tolerance = 0.1 + }; + + const auto lattices = Refine(spots, spots.size(), cells, scores, parameters); + CHECK(lattices.empty()); + for (const auto &l: lattices) + CHECK_NOTHROW(l.GetUBMatrix()); +}