From 02ecf5c32e93fb803495d0ae916acdb3e1e964fd Mon Sep 17 00:00:00 2001 From: leonarski_f Date: Fri, 5 Jun 2026 11:36:23 +0200 Subject: [PATCH] MultiLatticeSearch: Explore all valid sign flips when comparing with reference cell --- common/CrystalLattice.cpp | 9 ++- common/CrystalLattice.h | 1 + .../indexing/MultiLatticeSearch.cpp | 68 +++++++++++++------ tests/MultiLatticeSearchTest.cpp | 15 ++++ 4 files changed, 73 insertions(+), 20 deletions(-) diff --git a/common/CrystalLattice.cpp b/common/CrystalLattice.cpp index 0172b8fa..17dbca65 100644 --- a/common/CrystalLattice.cpp +++ b/common/CrystalLattice.cpp @@ -80,9 +80,16 @@ void CrystalLattice::Sort() { std::swap(vec[0], vec[1]); } +void CrystalLattice::FlipSign(size_t i1) { + if (i1 >= 3) + throw JFJochException(JFJochExceptionCategory::InputParameterInvalid, + "index out of range (0..2)"); + vec[i1] *= -1; +} + void CrystalLattice::FixHandedness() { if (CalcVolume() < 0) - vec[2] *= -1; + FlipSign(2); } diff --git a/common/CrystalLattice.h b/common/CrystalLattice.h index 5f4e6519..cb18b179 100644 --- a/common/CrystalLattice.h +++ b/common/CrystalLattice.h @@ -39,6 +39,7 @@ public: void Regularize(const gemmi::CrystalSystem &input); [[nodiscard]] std::vector GetUBMatrix() const; CrystalLattice NiggliReduce() const; + void FlipSign(size_t i1); }; inline std::ostream &operator<<( std::ostream &output, const CrystalLattice &in ) { diff --git a/image_analysis/indexing/MultiLatticeSearch.cpp b/image_analysis/indexing/MultiLatticeSearch.cpp index c799d9ab..11d0cf55 100644 --- a/image_analysis/indexing/MultiLatticeSearch.cpp +++ b/image_analysis/indexing/MultiLatticeSearch.cpp @@ -35,6 +35,30 @@ namespace { } } +CrystalLattice Transform(const CrystalLattice &in_latt, int i) { + CrystalLattice latt = in_latt; + switch (i) { + case 0: + break; + case 1: + latt.FlipSign(0); + latt.FlipSign(1); + break; + case 2: + latt.FlipSign(0); + latt.FlipSign(2); + break; + case 3: + latt.FlipSign(1); + latt.FlipSign(2); + break; + default: + break; + } + + return latt; +} + std::vector MultiLatticeSearch(const std::vector &lattices, float dist_tolerance, float angle_tolerance_deg) { @@ -48,29 +72,35 @@ std::vector MultiLatticeSearch(const std::vector(rod.x()), - static_cast(rod.y()), - static_cast(rod.z())); + const Eigen::AngleAxisd aa(R); + const Eigen::Vector3d rod = aa.angle() * aa.axis(); + const Coord rotation_vector(static_cast(rod.x()), + static_cast(rod.y()), + static_cast(rod.z())); - // output_lattice = R * reference, built straight from the Eigen matrix - // so it is exactly a proper rotation of the reference (full double precision). - const Eigen::Matrix3d out = R * LatticeMatrix(reference); - ret.push_back({ - latt, - CrystalLattice(Coord(out(0, 0), out(1, 0), out(2, 0)), - Coord(out(0, 1), out(1, 1), out(2, 1)), - Coord(out(0, 2), out(1, 2), out(2, 2))), - rotation_vector - }); + // output_lattice = R * reference, built straight from the Eigen matrix + // so it is exactly a proper rotation of the reference (full double precision). + const Eigen::Matrix3d out = R * LatticeMatrix(reference); + ret.push_back({ + latt, + CrystalLattice(Coord(out(0, 0), out(1, 0), out(2, 0)), + Coord(out(0, 1), out(1, 1), out(2, 1)), + Coord(out(0, 2), out(1, 2), out(2, 2))), + rotation_vector + }); + found = true; + } } return ret; diff --git a/tests/MultiLatticeSearchTest.cpp b/tests/MultiLatticeSearchTest.cpp index b40cc537..6e5df735 100644 --- a/tests/MultiLatticeSearchTest.cpp +++ b/tests/MultiLatticeSearchTest.cpp @@ -65,4 +65,19 @@ TEST_CASE("MultiLatticeSearch_SkipsDifferentCell") { TEST_CASE("MultiLatticeSearch_Empty") { auto result = MultiLatticeSearch({}); CHECK(result.empty()); +} + +TEST_CASE("MultiLatticeSearch_EP") { + // Real EP case + CrystalLattice cell1(Coord(-13.2, -30.0, -29.6), + Coord(70.1, -12.16, -18.6), + Coord(7.6, -23.9, 44.8)); + CrystalLattice cell2(Coord(-13.2, -29.9, -29.5), + Coord(-70.1, 12.15, 18.8), + Coord(1.8, 45.4, -23.7) + ); + + auto result = MultiLatticeSearch({cell1, cell2}, + 0.1, 5); + REQUIRE(result.size() == 2); } \ No newline at end of file