rotation_indexer: map FFT primitive cell to the SG conventional setting

When a space group is supplied without a reference cell, de-novo two-pass
rotation indexing fed the FFT's Niggli-reduced primitive cell straight into
XtalOptimizer as if it were the conventional cell. For non-primitive lattices
(centered I/F/R/C, or hexagonal where the primitive pair sits at gamma=60) the
conventional-system model then refined to a wrong minimum and indexed 0% of
frames: cytC (P3121) gave 103.9/103.9/78 instead of 83.7/83.7/88.6, insulin
(I213) 66.7 instead of 77.65, insulin-R3 51/51/36 instead of 81.4/81.4/33.3.

Run LatticeSearch on the FFT primitive cell (it already yields the correct
conventional cell + reindex for I/F/R/C). For the one remaining gap - a
metrically hexagonal lattice that the geometry-keyed search lands on the
ortho-hexagonal C setting - re-express the reduced primitive cell in
conventional hexagonal axes (b -> b - a opens gamma 60 -> 120).

De-novo "-S" now indexes cytC/insu/Ins_H/lyso/MyoB/EP/lyso_ref at 100% with the
correct cell; the "-C -S" path is unchanged. The helper stays in this .cpp
(g++) rather than CrystalLattice.h to avoid recompiling CUDA units, which is
broken under the box's CUDA-13 nvcc; promote it to a method once that is fixed.

Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
This commit is contained in:
2026-06-30 14:29:50 +02:00
co-authored by Claude Opus 4.8
parent 3aec235f67
commit 21dd99524f
@@ -9,6 +9,19 @@
#include "../lattice_search/LatticeSearch.h"
#include "../indexing/MultiLatticeSearch.h"
namespace {
// Re-express a primitive hexagonal/trigonal lattice in the conventional hexagonal setting
// (a = b, gamma = 120). The Niggli-reduced primitive cell carries the two equal-length axes
// at gamma = 60; replacing b with b - a opens that angle to 120 without changing the lattice.
CrystalLattice HexagonalConventional(CrystalLattice latt) {
latt.ReorderABEqual(); // put the equal-length pair in a, b
Coord a = latt.Vec0(), b = latt.Vec1(), c = latt.Vec2();
if (angle_deg(a, b) < 90.0f)
b -= a;
return CrystalLattice(a, b, c); // constructor fixes handedness
}
}
RotationIndexer::RotationIndexer(const DiffractionExperiment &x, IndexerThreadPool &indexer)
: experiment(x),
index_ice_rings(x.GetIndexingSettings().GetIndexIceRings()),
@@ -40,12 +53,28 @@ void RotationIndexer::RunIndexing() {
if (!indexer_result.lattice.empty() && indexer_result.lattice[0].CalcVolume() > 1.0) {
auto sg = experiment.GetGemmiSpaceGroup();
if (sg) {
// The FFT returns a primitive cell; map it to the SG's conventional setting via the
// Niggli-class table (which handles I/F/R/C centering). LatticeSearch keys off geometry
// alone, so a metrically-hexagonal lattice can land on the ortho-hexagonal C setting -
// when the SG is trigonal/hexagonal, re-express the reduced primitive cell in the
// conventional hexagonal axes instead.
auto ls = LatticeSearch(indexer_result.lattice[0]);
const auto is_hexagonal = [](gemmi::CrystalSystem s) {
return s == gemmi::CrystalSystem::Trigonal || s == gemmi::CrystalSystem::Hexagonal;
};
CrystalLattice conventional = ls.conventional;
if (is_hexagonal(sg->crystal_system()) && !is_hexagonal(ls.system))
conventional = HexagonalConventional(ls.primitive_reduced);
search_result_ = LatticeSearchResult{
.niggli_class = 0, // Since Niggli class was not searched for, we don't know which one
.conventional = indexer_result.lattice[0], // If lattice provided, it is for now primitive == conventional
.niggli_class = ls.niggli_class,
.primitive_reduced = ls.primitive_reduced,
.conventional = conventional,
.system = sg->crystal_system(),
.centering = sg->centring_type(),
.reindex = gemmi::Mat33(1, 0, 0, 0, 1, 0, 0, 0, 1),
.reindex = ls.reindex,
};
} else {
// Find lattice type based on cell