diff --git a/image_analysis/rotation_indexer/RotationIndexer.cpp b/image_analysis/rotation_indexer/RotationIndexer.cpp index f01c8f60..1276ba2b 100644 --- a/image_analysis/rotation_indexer/RotationIndexer.cpp +++ b/image_analysis/rotation_indexer/RotationIndexer.cpp @@ -274,6 +274,7 @@ void RotationIndexer::RunIndexing() { size_t best_ci = 0; XtalOptimizerData best_data; LatticeSearchResult best_sr; + std::shared_ptr best_alt; for (size_t ci = 0; ci < n_try; ci++) { if (!work[ci].viable) continue; @@ -290,16 +291,33 @@ void RotationIndexer::RunIndexing() { // symmetry (incl. R-centred) indexes comparably (ratio ~0.7). Preferring the constrained // cell on a near-tie keeps the real symmetry/centering; the intensities settle the final // space group. + // When the guard leaves the constrained cell in place, hand the refined triclinic cell to + // the caller instead of dropping it: the accumulated-spot fraction separates a false + // promotion from genuine symmetry by less than 2x on a lattice that is pseudo-symmetric + // to a few tenths of a degree, while the caller's per-frame validation separates the same + // pair by more than 20x. + std::shared_ptr tri_alt; if (work[ci].has_tri) { Solved t = tri_f[ci].get(); + auto as_triclinic = [](LatticeSearchResult s) { + s.system = gemmi::CrystalSystem::Triclinic; + s.centering = 'P'; + s.conventional = s.primitive_reduced; + s.reindex = gemmi::Mat33(1, 0, 0, 0, 1, 0, 0, 0, 1); + return s; + }; if (t.ok && t.frac > 0.3f && frac < 0.5f * t.frac) { + sr = as_triclinic(sr); data = std::move(t.data); ok = true; frac = t.frac; - sr.system = gemmi::CrystalSystem::Triclinic; - sr.centering = 'P'; - sr.conventional = sr.primitive_reduced; - sr.reindex = gemmi::Mat33(1, 0, 0, 0, 1, 0, 0, 0, 1); + } else if (t.ok) { + tri_alt = std::make_shared(RotationIndexerResult{ + .lattice = t.data.latt, + .search_result = as_triclinic(sr), + .geom = t.data.geom, + .axis = t.data.axis, + }); } } @@ -364,6 +382,7 @@ void RotationIndexer::RunIndexing() { best_data = std::move(data); best_sr = sr; best_ci = ci; + best_alt = std::move(tri_alt); } } @@ -372,6 +391,7 @@ void RotationIndexer::RunIndexing() { indexed_lattice = best_data.latt; updated_geom_ = best_data.geom; axis_ = best_data.axis; + unconstrained_ = std::move(best_alt); } // Extra (twin) lattices: MultiLatticeSearch derives each rotation by relating the FFT's primary @@ -469,6 +489,7 @@ std::optional RotationIndexer::GetLattice() const { .search_result = search_result_, .geom = updated_geom_, .axis = axis_, + .unconstrained = unconstrained_, }; } @@ -479,6 +500,7 @@ void RotationIndexer::ForceResult(const RotationIndexerResult &result) { search_result_ = result.search_result; updated_geom_ = result.geom; axis_ = result.axis; + unconstrained_ = result.unconstrained; } bool RotationIndexer::AccumulationFull() const { diff --git a/image_analysis/rotation_indexer/RotationIndexer.h b/image_analysis/rotation_indexer/RotationIndexer.h index c0049017..1fe5de58 100644 --- a/image_analysis/rotation_indexer/RotationIndexer.h +++ b/image_analysis/rotation_indexer/RotationIndexer.h @@ -3,6 +3,7 @@ #pragma once +#include #include #include @@ -17,6 +18,12 @@ struct RotationIndexerResult { LatticeSearchResult search_result; DiffractionGeometry geom; std::optional axis; + // The same lattice refined WITHOUT the metric symmetry (triclinic P), when there is one to + // compare against. The Bravais class is decided on the UNREFINED FFT candidate against a fixed + // angular tolerance, so a pseudo-symmetric cell can be promoted a class too far and the + // constraint then snaps a real angle to the ideal one; here the caller can settle that on a + // statistic the accumulated-spot fraction is too blunt for. Null on the alternative itself. + std::shared_ptr unconstrained; }; class RotationIndexer { @@ -38,6 +45,7 @@ class RotationIndexer { DiffractionGeometry updated_geom_; LatticeSearchResult search_result_; std::vector extra_lattices_; + std::shared_ptr unconstrained_; IndexerThreadPool &indexer_; diff --git a/rugnux/Rugnux.cpp b/rugnux/Rugnux.cpp index 8ef1e368..f069209a 100644 --- a/rugnux/Rugnux.cpp +++ b/rugnux/Rugnux.cpp @@ -710,10 +710,47 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b for (size_t i = 0; i < ris.size(); i++) { if (cancelled_) break; - const auto result = ris[i]->GetLattice(); - if (!result.has_value()) + auto found = ris[i]->GetLattice(); + if (!found.has_value()) continue; - const int score = count_indexed(idx, *result); + RotationIndexerResult result = std::move(*found); + int score = count_indexed(idx, result); + const std::string &name = schemes[i].first; + logger.Info("First-pass scheme '{}': indexes {}/{} validation frames", name, score, + static_cast(validation.size())); + + // Metric symmetry that costs nearly every frame is not this crystal's symmetry. The + // Bravais class is decided on the UNREFINED FFT candidate against a fixed angular + // tolerance, so a lattice that is pseudo-symmetric to a few tenths of a degree gets + // promoted a class too far; the constraint then snaps a real angle to the ideal one, + // which throws nearly every reflection of every frame out of tolerance (measured: 2/60 + // frames for a monoclinic lattice promoted to C-orthorhombic, 39/60 for its primitive + // cell). The indexer's own guard compares the two on the accumulated first-pass cloud, + // where they differ by less than a factor 2 - far too little to act on. + // The bar is a clear majority, not a margin: the unconstrained refinement holds NO cell + // parameter fixed, so it can only index at least as many frames as the constrained one, + // and on genuine symmetry it does index a few more (measured: 47->54 frames on a real + // I-centred orthorhombic crystal, 49->58 on an F-cubic one - a 10% margin demotes both). + // Only a constrained cell that fails outright while its unconstrained cell works is + // evidence of a false promotion, so demand exactly that. It is the same "fails to index + // half the frames" test the long-axis rescue below uses, and it leaves genuine symmetry + // alone; the space group is in any case decided from the intensities later. + const int majority = static_cast(validation.size()) / 2; + if (result.unconstrained && score < majority) { + RotationIndexerResult alt = *result.unconstrained; + const int alt_score = count_indexed(idx, alt); + if (alt_score > majority) { + logger.Info("Scheme '{}': {}-centred {} indexes {}/{} frames but its unconstrained " + "cell indexes {}/{} - the metric symmetry is a false promotion, dropping it", + name, result.search_result.centering, + gemmi::crystal_system_str(result.search_result.system), score, + static_cast(validation.size()), alt_score, + static_cast(validation.size())); + result = std::move(alt); + score = alt_score; + } + } + // Compare PRIMITIVE volumes: two schemes can find the same lattice in different // settings, and a centred setting's cell is an exact integer multiple of its primitive // one - a rhombohedral lattice in hexagonal axes is exactly 3x its primitive @@ -721,10 +758,7 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b // on that pair and "demote" a perfectly good setting to a threefold-smaller merge, which // is enough to change the space group the search then picks. const double vol = std::abs( - result->lattice.ToPrimitive(result->search_result.centering).CalcVolume()); - const std::string &name = schemes[i].first; - logger.Info("First-pass scheme '{}': indexes {}/{} validation frames", name, score, - static_cast(validation.size())); + result.lattice.ToPrimitive(result.search_result.centering).CalcVolume()); // A later scheme wins if it indexes clearly more frames (>10%). const bool clearly_more = static_cast(score) > bp.score * 1.1f + 0.5f; @@ -754,7 +788,7 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b bp.score = score; bp.vol = vol; bp.name = name; - bp.result = result; + bp.result = std::move(result); } } return bp;