diff --git a/image_analysis/rotation_indexer/RotationIndexer.cpp b/image_analysis/rotation_indexer/RotationIndexer.cpp index 2b40f7e7..4334abb9 100644 --- a/image_analysis/rotation_indexer/RotationIndexer.cpp +++ b/image_analysis/rotation_indexer/RotationIndexer.cpp @@ -20,6 +20,27 @@ namespace { b -= a; return CrystalLattice(a, b, c); // constructor fixes handedness } + + bool IsHexagonalSystem(gemmi::CrystalSystem s) { + return s == gemmi::CrystalSystem::Trigonal || s == gemmi::CrystalSystem::Hexagonal; + } + + // The hexagonal lattice metric (two equal axes at 60/120 deg, both perpendicular to the third) is + // also satisfied by its ortho-hexagonal C-centred supercell, so the geometry-keyed LatticeSearch can + // land there. Detect the hexagonal metric on the reduced PRIMITIVE cell so the de-novo path (no space + // group to key on) can re-express it in conventional hexagonal axes. + bool IsMetricallyHexagonal(CrystalLattice latt, float rel_tol = 0.03f, float angle_tol_deg = 3.0f) { + latt.ReorderABEqual(); + const Coord a = latt.Vec0(), b = latt.Vec1(), c = latt.Vec2(); + const float la = a.Length(), lb = b.Length(); + if (la <= 0.0f || lb <= 0.0f || std::fabs(la - lb) > rel_tol * std::max(la, lb)) + return false; + const float gab = angle_deg(a, b); + if (std::fabs(gab - 60.0f) > angle_tol_deg && std::fabs(gab - 120.0f) > angle_tol_deg) + return false; + return std::fabs(angle_deg(a, c) - 90.0f) <= angle_tol_deg && + std::fabs(angle_deg(b, c) - 90.0f) <= angle_tol_deg; + } } RotationIndexer::RotationIndexer(const DiffractionExperiment &x, IndexerThreadPool &indexer) @@ -77,8 +98,16 @@ void RotationIndexer::RunIndexing() { .reindex = ls.reindex, }; } else { - // Find lattice type based on cell - search_result_ = LatticeSearch(indexer_result.lattice[0]); + // Find lattice type from the metric. LatticeSearch keys off geometry, so a metrically + // hexagonal lattice can land on the ortho-hexagonal C setting; re-express it in conventional + // hexagonal axes (P, gamma=120) so the trigonal/hexagonal 3-fold is not hidden from scaling. + auto ls = LatticeSearch(indexer_result.lattice[0]); + if (!IsHexagonalSystem(ls.system) && IsMetricallyHexagonal(ls.primitive_reduced)) { + ls.conventional = HexagonalConventional(ls.primitive_reduced); + ls.system = gemmi::CrystalSystem::Hexagonal; + ls.centering = 'P'; + } + search_result_ = ls; } // Run refinement