diff --git a/docs/CHANGELOG.md b/docs/CHANGELOG.md index 306b4aab..9f77680e 100644 --- a/docs/CHANGELOG.md +++ b/docs/CHANGELOG.md @@ -22,6 +22,8 @@ This is an UNSTABLE release. It includes many experimental features, as well as * rugnux: scaling and merging are faster, with identical output. * The per-image resolution estimate now predicts the resolution the merged data reach, rather than reporting the highest-resolution spot found; rugnux reports the run's value as `SPOT_RESOLUTION_ESTIMATE` in its report. * rugnux: fixing the space group with `-S` no longer prevents the lattice from being found; the group is applied to scaling and merging rather than to the indexing search. +* rugnux: a lattice indexed in a different setting from the space group fixed with `-S` is reindexed into that group's own setting before merging, instead of merging with the group's centring rule applied in the wrong frame. +* rugnux: a run now stops, naming the cell the crystal indexed as, when the space group fixed with `-S` has a lattice the crystal does not have, rather than merging and reporting statistics that cannot describe it. * rugnux: the detector geometry is also logged in XDS's convention (`ORGX`/`ORGY`, detector axis vectors, rotation axis), so it can be compared with an XDS refinement. * rugnux: `_process.h5` describes the pixel format of the images it links to. * A DECTRIS detector sending signed images is no longer declared unsigned in the image stream and in HDF5. diff --git a/docs/RUGNUX.md b/docs/RUGNUX.md index a1b910c4..c7aac77a 100644 --- a/docs/RUGNUX.md +++ b/docs/RUGNUX.md @@ -328,7 +328,9 @@ and both are written into the merged `.cif`. Run **fully de novo** (no `-C`/`-S`) for the best result — supplying a cell or space group up front -tends to *degrade* low-symmetry cases. `--scaling-high-resolution` (set it to your expected +tends to *degrade* low-symmetry cases. A `-S` group whose Bravais lattice the crystal turns out not to +have stops the run and names the cell that was indexed, rather than merging in a frame the reflections +are not in; where the lattice does have that group's setting, the reflections are reindexed into it. `--scaling-high-resolution` (set it to your expected resolution) sharpens both the space-group search and the error model. To tune the first pass use `--two-pass-rotation=100` (or `-R100` — the first-pass image count); to force the sweep to be treated as independent stills use `--force-still`. diff --git a/image_analysis/IndexAndRefine.cpp b/image_analysis/IndexAndRefine.cpp index 7e43d049..c3df55bd 100644 --- a/image_analysis/IndexAndRefine.cpp +++ b/image_analysis/IndexAndRefine.cpp @@ -499,8 +499,11 @@ void IndexAndRefine::QuickPredictAndIntegrate(DataMessage &msg, // Centering is a hypothesis to confirm, not assume: with no user-fixed space group, predict // in P so the centering-absent reflections are integrated and the space-group search can // confirm or disprove centering (and catch a missed superstructure). A user-fixed space - // group is trusted, so reject its absences here. - .centering = experiment.GetGemmiSpaceGroup().has_value() ? outcome.symmetry.centering : 'P', + // group is trusted, so reject its absences here - unless this pass measures geometry and + // discards its intensities (predict_all_centring_nodes_), where the absences only halve the + // events the geometry is fitted from. + .centering = experiment.GetGemmiSpaceGroup().has_value() && !predict_all_centring_nodes_ + ? outcome.symmetry.centering : 'P', .wedge_deg = std::fabs(wedge_deg), .mosaicity_deg = std::fabs(mos_deg), // FWHM -> sigma; 0 when monochromatic, leaving the prediction unchanged. diff --git a/image_analysis/IndexAndRefine.h b/image_analysis/IndexAndRefine.h index d10e66c7..fb3f955e 100644 --- a/image_analysis/IndexAndRefine.h +++ b/image_analysis/IndexAndRefine.h @@ -69,6 +69,13 @@ class IndexAndRefine { // uses the frame-order-SMOOTHED mosaicity that RotationScaleMerge already fitted in the first pass, // rather than re-deriving it from scratch. std::vector prediction_mosaicity_override_; + // Predict every node of the lattice, ignoring the centring absences of a fixed space group. Set + // for the rotation two-pass GEOMETRY pre-pass, whose job is to measure the detector geometry from + // spot positions and whose intensities are thrown away: rejecting the absences there costs it half + // its events and buys nothing. Measured with an I-centred group fixed - the pre-pass fitted a + // different error model (ISa 7.8 -> 3.6), post-refined the distance 119 um away, and the second + // pass re-indexed 49 of 60 frames instead of 60. + bool predict_all_centring_nodes_ = false; std::vector scale_cc; std::vector > unit_cells; @@ -106,6 +113,8 @@ public: void SetPredictionMosaicityOverride(std::vector mosaicity_per_frame) { prediction_mosaicity_override_ = std::move(mosaicity_per_frame); } + // Predict the centring-absent reflections too, even with a fixed space group - see the member. + void PredictAllCentringNodes(bool on) { predict_all_centring_nodes_ = on; } // Returns whether the frame indexed (a lattice was found and refined). Integration, when it runs, // is a further step gated on quick_integration. diff --git a/image_analysis/lattice_search/LatticeSearch.cpp b/image_analysis/lattice_search/LatticeSearch.cpp index d8e5e85a..ed6e7df7 100644 --- a/image_analysis/lattice_search/LatticeSearch.cpp +++ b/image_analysis/lattice_search/LatticeSearch.cpp @@ -32,7 +32,13 @@ struct NiggliClass { char centering; }; -LatticeSearchResult LatticeSearch(const CrystalLattice &L, double dist_tolerance, double angle_tolerance) { +namespace { +// The body of both entry points. only_class, when given, keeps just the characters of that Bravais +// class - see LatticeSearchForClass. With no filter this is the original walk unchanged, and the +// triclinic character fits every metric, so it always returns a result. +std::optional SearchCharacters(const CrystalLattice &L, double dist_tolerance, + double angle_tolerance, + const std::pair *only_class) { UnitCell uc = L.GetUnitCell(); gemmi::UnitCell g_uc(uc.a, uc.b, uc.c, uc.alpha, uc.beta, uc.gamma); gemmi::GruberVector g_vec(g_uc, 'P', true); @@ -320,6 +326,8 @@ LatticeSearchResult LatticeSearch(const CrystalLattice &L, double dist_tolerance -> std::optional { const auto uc_reduced = latt.GetUnitCell(); for (const auto &c: make_classes(D, E, F)) { + if (only_class && (c.system != only_class->first || c.centering != only_class->second)) + continue; if (c.type == 1 && uc_reduced.beta >= 90 - angle_tolerance ) continue; @@ -384,6 +392,9 @@ LatticeSearchResult LatticeSearch(const CrystalLattice &L, double dist_tolerance if (found) return *found; + if (only_class) + return std::nullopt; // no character of the requested class fits this metric + return LatticeSearchResult{ .niggli_class = 44, .primitive_reduced = L_niggli, @@ -393,3 +404,15 @@ LatticeSearchResult LatticeSearch(const CrystalLattice &L, double dist_tolerance .reindex = gemmi::Mat33(1, 0, 0, 0, 1, 0, 0, 0, 1), }; } +} // namespace + +LatticeSearchResult LatticeSearch(const CrystalLattice &L, double dist_tolerance, double angle_tolerance) { + return *SearchCharacters(L, dist_tolerance, angle_tolerance, nullptr); +} + +std::optional LatticeSearchForClass(const CrystalLattice &L, + gemmi::CrystalSystem system, char centering, + double dist_tolerance, double angle_tolerance) { + const std::pair only{system, centering}; + return SearchCharacters(L, dist_tolerance, angle_tolerance, &only); +} diff --git a/image_analysis/lattice_search/LatticeSearch.h b/image_analysis/lattice_search/LatticeSearch.h index 422d1c50..351dd762 100644 --- a/image_analysis/lattice_search/LatticeSearch.h +++ b/image_analysis/lattice_search/LatticeSearch.h @@ -7,6 +7,7 @@ #include "../../common/UnitCell.h" #include "../../common/Coord.h" #include +#include #include "gemmi/symmetry.hpp" struct LatticeSearchResult { @@ -23,3 +24,15 @@ struct LatticeSearchResult { LatticeSearchResult LatticeSearch(const CrystalLattice& L, double dist_tolerance = 0.03, double angle_tolerance_deg = 3); +// The same 44-character walk, restricted to the characters of ONE Bravais class: the setting of L +// that carries exactly (system, centering), or nothing when no character of that class fits the +// metric. LatticeSearch returns the FIRST character that fits, which is the most symmetric setting +// the metric supports, so a tetragonal-P lattice never comes back in its C-centred orthorhombic +// setting even though the same lattice has one. Used where a space group has been fixed by hand and +// the setting the reflections must be indexed in is that group's own lattice, not the most +// symmetric one. A trigonal-P group has a hexagonal-P lattice - ask for Hexagonal, 'P'. +std::optional LatticeSearchForClass(const CrystalLattice& L, + gemmi::CrystalSystem system, char centering, + double dist_tolerance = 0.03, + double angle_tolerance_deg = 3); + diff --git a/rugnux/Rugnux.cpp b/rugnux/Rugnux.cpp index 57bbd2ff..6b10ffe3 100644 --- a/rugnux/Rugnux.cpp +++ b/rugnux/Rugnux.cpp @@ -266,7 +266,8 @@ namespace { Rugnux::Rugnux(JFJochHDF5Reader &reader, DiffractionExperiment experiment, PixelMask pixel_mask, ProcessConfig config) : reader_(reader), experiment_(std::move(experiment)), - pixel_mask_(std::move(pixel_mask)), config_(std::move(config)) { + pixel_mask_(std::move(pixel_mask)), config_(std::move(config)), + user_fixed_sg_(experiment_.GetSpaceGroupNumber()) { // Bit 9 describes where THIS run found the beam stop, so a mask read back from a file that // already carries one starts clear; the user mask (bit 8) is left as it was loaded. pixel_mask_.ClearBeamStopMask(experiment_); @@ -1235,6 +1236,11 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b // Second pass of the rotation two-pass: predict with the smoothed mosaicity fitted in the pre-pass. if (!geometry_prepass && !prepass_mosaicity_.empty()) indexer->SetPredictionMosaicityOverride(prepass_mosaicity_); + // The geometry pre-pass fits the detector geometry to spot POSITIONS and throws its intensities + // away, so a fixed group's centring absences cost it half its events for nothing. A no-op when + // no group is fixed - prediction is already in P there - so the de-novo path is untouched. + if (geometry_prepass) + indexer->PredictAllCentringNodes(true); } const auto start_time = std::chrono::steady_clock::now(); @@ -2034,6 +2040,43 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b // (forced) mosaicity is handled by the recompute. The decay + absorption correction surfaces run as // post-scale-fulls stages (on by default, see ScalingSettings::CorrectionSurfaces); external- // reference scaling, the stills -B (per-image B-factor) and wedge refinement are unsupported here. + // Re-seat the integrated reflections into a different setting of the SAME lattice, for the two + // arms below. The change of basis is read straight off the two lattices: P[i][j] = + // conv_real[i] . indexed_reciprocal[j] is the coefficient of the indexed cell's j-th axis in the + // new cell's i-th axis, and indices transform with the axes, so it is also the matrix that takes + // hkl across. Rounding it to integers - and declining when the residual is large - keeps it exact + // even though the refined cell is not exactly LatticeSearch's Niggli cell. + const auto reindex_into = [&](const LatticeSearchResult &cand) { + const Coord cv[3] = {cand.conventional.Vec0(), cand.conventional.Vec1(), cand.conventional.Vec2()}; + const Coord rs[3] = {end_msg.rotation_lattice->Astar(), end_msg.rotation_lattice->Bstar(), + end_msg.rotation_lattice->Cstar()}; + gemmi::Mat33 reindex; + double reindex_res = 0.0; + for (int i = 0; i < 3; ++i) + for (int j = 0; j < 3; ++j) { + const double v = cv[i] * rs[j]; + reindex.a[i][j] = std::round(v); + reindex_res = std::max(reindex_res, std::fabs(v - reindex.a[i][j])); + } + if (reindex_res >= 0.1) + return false; + for (auto &io : indexer->GetIntegrationOutcome()) { + io.latt = io.latt.Multiply(reindex); + for (auto &r : io.reflections) { + const double h = r.h, k = r.k, l = r.l; + r.h = static_cast(std::lround(reindex.a[0][0]*h + reindex.a[0][1]*k + reindex.a[0][2]*l)); + r.k = static_cast(std::lround(reindex.a[1][0]*h + reindex.a[1][1]*k + reindex.a[1][2]*l)); + r.l = static_cast(std::lround(reindex.a[2][0]*h + reindex.a[2][1]*k + reindex.a[2][2]*l)); + } + } + result.consensus_cell = cand.conventional.GetUnitCell(); + end_msg.unit_cell = result.consensus_cell; + end_msg.rotation_lattice = end_msg.rotation_lattice->Multiply(reindex); + end_msg.rotation_lattice_type = LatticeMessage{ .centering = cand.centering, + .niggli_class = end_msg.rotation_lattice_type->niggli_class, .crystal_system = cand.system }; + return true; + }; + // Two-pass second pass: if the de-novo indexer DEMOTED to a triclinic (primitive) cell while the reused // merge space group is higher-symmetry (a huge oblique cell whose constrained centred refine was // ill-posed, so the pseudo-symmetry guard kept the primitive), the reflections are in the primitive @@ -2045,39 +2088,44 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b const auto *msg_sg = gemmi::find_spacegroup_by_number(static_cast(*prepass_merge_sg_)); if (msg_sg && msg_sg->crystal_system() != gemmi::CrystalSystem::Triclinic) { const auto cand = LatticeSearch(*end_msg.rotation_lattice); - // Exact integer reindex prim->conv: P[i][j] = conv_real[i] . rot_reciprocal[j] (as the - // centred-lattice test), robust to the refined cell not being exactly LatticeSearch's Niggli cell. - const Coord cv[3] = {cand.conventional.Vec0(), cand.conventional.Vec1(), cand.conventional.Vec2()}; - const Coord rs[3] = {end_msg.rotation_lattice->Astar(), end_msg.rotation_lattice->Bstar(), - end_msg.rotation_lattice->Cstar()}; - gemmi::Mat33 reindex; - double reindex_res = 0.0; - for (int i = 0; i < 3; ++i) - for (int j = 0; j < 3; ++j) { - const double v = cv[i] * rs[j]; - reindex.a[i][j] = std::round(v); - reindex_res = std::max(reindex_res, std::fabs(v - reindex.a[i][j])); - } // Only when LatticeSearch recovers the reused group's own metric (system + centring): otherwise // the reused group and the indexed metric disagree and reindexing would be wrong - leave it. if (cand.system == msg_sg->crystal_system() && cand.centering == msg_sg->centring_type() - && reindex_res < 0.1) { - for (auto &io : indexer->GetIntegrationOutcome()) { - io.latt = io.latt.Multiply(reindex); - for (auto &r : io.reflections) { - const double h = r.h, k = r.k, l = r.l; - r.h = static_cast(std::lround(reindex.a[0][0]*h + reindex.a[0][1]*k + reindex.a[0][2]*l)); - r.k = static_cast(std::lround(reindex.a[1][0]*h + reindex.a[1][1]*k + reindex.a[1][2]*l)); - r.l = static_cast(std::lround(reindex.a[2][0]*h + reindex.a[2][1]*k + reindex.a[2][2]*l)); - } - } - result.consensus_cell = cand.conventional.GetUnitCell(); - end_msg.unit_cell = result.consensus_cell; - end_msg.rotation_lattice = end_msg.rotation_lattice->Multiply(reindex); - end_msg.rotation_lattice_type = LatticeMessage{ .centering = cand.centering, - .niggli_class = end_msg.rotation_lattice_type->niggli_class, .crystal_system = cand.system }; + && reindex_into(cand)) logger.Info("Two-pass: reindexed the de-novo primitive cell into the space-group {} " "conventional setting for the merge", static_cast(*prepass_merge_sg_)); + } + } + + // A space group the USER fixed carries its own Bravais lattice, and that is the setting its + // reflections have to be indexed in - but the indexer answers to the metric, not to the group, and + // LatticeSearch hands back the MOST symmetric setting the metric supports. A tetragonal-P lattice + // therefore never comes back C-centred orthorhombic even though the same lattice has that setting, + // and -S C222 on one merges with the C absence rule applied in the wrong frame: measured, half the + // observations thrown away (2.14 M -> 1.08 M) and ISa 28.3 -> 24.6. Ask the character table for the + // group's own class instead and re-seat the reflections into it. + // + // Keyed on the CENTRING differing, which is the whole of the harm: a fixed group of lower symmetry + // than the lattice but with the same centring (-S P21 on an orthorhombic-P lattice, say) merges + // correctly where it stands, and re-seating it would only permute axes that are already right. + // A group the RUN determined is left alone - the arm above and the centring check below are its + // handling, and a de-novo run must behave exactly as it did. + if (user_fixed_sg_ && end_msg.rotation_lattice.has_value() && end_msg.rotation_lattice_type.has_value()) { + const auto *fixed_sg = gemmi::find_spacegroup_by_number(static_cast(*user_fixed_sg_)); + const char indexed_centering = end_msg.rotation_lattice_type->centering; + if (fixed_sg && indexed_centering != fixed_sg->centring_type()) { + // A trigonal-P group sits on a hexagonal-P lattice; every other system names its own. + const auto want_system = fixed_sg->crystal_system() == gemmi::CrystalSystem::Trigonal + && fixed_sg->centring_type() == 'P' + ? gemmi::CrystalSystem::Hexagonal : fixed_sg->crystal_system(); + const auto cand = LatticeSearchForClass(*end_msg.rotation_lattice, want_system, + fixed_sg->centring_type()); + if (cand && reindex_into(*cand)) { + const auto &uc = *result.consensus_cell; + logger.Info("Reindexed the {}-centred indexed lattice into the {}-centred setting the " + "fixed space group {} needs: a={:.3f} b={:.3f} c={:.3f} alpha={:.2f} " + "beta={:.2f} gamma={:.2f}", indexed_centering, cand->centering, + fixed_sg->xhm(), uc.a, uc.b, uc.c, uc.alpha, uc.beta, uc.gamma); } } } @@ -2106,6 +2154,39 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b result.lattice_conflicts_with_prepass_sg = true; } } + + // The same conflict, but for a group the USER fixed, and on a pass that has nothing to fall back + // on - the first of the two, or a single pass. The re-seating above has had its chance and no + // setting of the lattice this crystal indexes as carries the group's Bravais lattice, so there is + // no frame in which the merge means anything. What comes out is not a merely suboptimal answer: + // the absence rule is applied across a frame the reflections are not in, and the statistics stop + // being arithmetic (measured: 173.5% complete on a triclinic-P lattice merged in C2, and an + // undefined R_meas on an F-centred cubic one merged in a trigonal-P group). Nor is there anything + // to fall back ON - the group is the user's assertion, and quietly determining a different one + // would answer a question that was not asked. Refuse, and name the cell that WAS indexed so the + // user can act on it. + // + // The SECOND pass keeps the flag-and-do-not-adopt handling just above instead: its lattice comes + // from a de-novo re-index at the post-refined geometry, and the first pass - whose lattice did + // carry the group - is still there to go back to. Refusing there would throw away a good answer. + if (user_fixed_sg_ && !prepass_merge_sg_ + && end_msg.rotation_lattice_type.has_value() && end_msg.unit_cell.has_value()) { + const auto *fixed_sg = gemmi::find_spacegroup_by_number(static_cast(*user_fixed_sg_)); + if (fixed_sg && end_msg.rotation_lattice_type->centering != fixed_sg->centring_type()) { + const auto &uc = *end_msg.unit_cell; + throw JFJochException(JFJochExceptionCategory::InputParameterInvalid, fmt::format( + "The space group {} was fixed for this run, and its lattice is {}-centred - but " + "this crystal indexes as a {}-centred {} lattice (a={:.3f} b={:.3f} c={:.3f} " + "alpha={:.2f} beta={:.2f} gamma={:.2f}), and no setting of that lattice carries the " + "group's. Merging in {} would apply its absence rule to reflections that are not in " + "its frame, so nothing it reported would describe this crystal. Re-run without a fixed " + "space group to have it determined from the data, or fix one whose lattice this " + "crystal has.", + fixed_sg->xhm(), fixed_sg->centring_type(), end_msg.rotation_lattice_type->centering, + gemmi::crystal_system_str(end_msg.rotation_lattice_type->crystal_system), + uc.a, uc.b, uc.c, uc.alpha, uc.beta, uc.gamma, fixed_sg->xhm())); + } + } if (prepass_merge_sg_) experiment_.SpaceGroupNumber(*prepass_merge_sg_); diff --git a/rugnux/Rugnux.h b/rugnux/Rugnux.h index e8fee939..2e6c40fe 100644 --- a/rugnux/Rugnux.h +++ b/rugnux/Rugnux.h @@ -269,6 +269,14 @@ class Rugnux { // the second pass, and for a reindex-derived group (left to re-search). std::optional prepass_merge_sg_; + // The space group the USER fixed (-S on the CLI, the settings panel in the viewer), captured before + // anything can clear it - the two-pass driver blanks the experiment's group between the passes so the + // second pass re-indexes de novo, so after that point the experiment can no longer say whether the + // group in force was asserted or determined. A group the user asserted is treated differently from + // one the run determined: it is re-seated onto its own Bravais lattice if the metric has that + // setting, and the run refuses rather than merging in it if the metric does not. + const std::optional user_fixed_sg_; + // Whether the group in prepass_merge_sg_ was PROMOTED by pass-1's search rather than given. Pass 2 // reuses the group without searching, so it cannot work this out for itself - and without it the // canonical output reports "the Laue class is holohedral, so no twin law exists" about a Laue class diff --git a/tests/LatticeSearchTest.cpp b/tests/LatticeSearchTest.cpp index 9957abec..0eb53bb3 100644 --- a/tests/LatticeSearchTest.cpp +++ b/tests/LatticeSearchTest.cpp @@ -7,6 +7,7 @@ #include "../common/UnitCell.h" #include "../image_analysis/lattice_search/LatticeSearch.h" #include "gemmi/symmetry.hpp" +#include // Helper: check near-equality of unit cell parameters static void check_uc(const UnitCell& uc, double a, double b, double c, @@ -459,3 +460,55 @@ TEST_CASE("LatticeSearch - trigonal R") { CHECK(uc_prim.beta == Catch::Approx(alpha).margin(1e-2)); CHECK(uc_prim.gamma == Catch::Approx(alpha).margin(1e-2)); } + +// The class-filtered walk: the same table, restricted to one Bravais class. A tetragonal-P lattice is +// also a C-centred orthorhombic one (a_C = a+b, b_C = -a+b, c_C = c), and asking for that class has to +// return that setting even though the plain search rightly prefers the tetragonal one. +TEST_CASE("LatticeSearchForClass - tetragonal P also has a C-centred orthorhombic setting") { + const double a = 50.0, c = 120.0; + const CrystalLattice L(a, a, c, 90, 90, 90); + + const auto plain = LatticeSearch(L, 1e-6); + CHECK(plain.system == gemmi::CrystalSystem::Tetragonal); + CHECK(plain.centering == 'P'); + + const auto ortho = LatticeSearchForClass(L, gemmi::CrystalSystem::Orthorhombic, 'C', 1e-6); + REQUIRE(ortho.has_value()); + CHECK(ortho->system == gemmi::CrystalSystem::Orthorhombic); + CHECK(ortho->centering == 'C'); + const auto uc = ortho->conventional.GetUnitCell(); + // The C cell is the face diagonal on a and b, so twice the volume and a = b = a_tet * sqrt(2). + CHECK(uc.a == Catch::Approx(a * std::sqrt(2.0)).margin(1e-4)); + CHECK(uc.b == Catch::Approx(a * std::sqrt(2.0)).margin(1e-4)); + CHECK(uc.c == Catch::Approx(c).margin(1e-4)); + CHECK(uc.alpha == Catch::Approx(90).margin(1e-4)); + CHECK(uc.beta == Catch::Approx(90).margin(1e-4)); + CHECK(uc.gamma == Catch::Approx(90).margin(1e-4)); +} + +TEST_CASE("LatticeSearchForClass - a class the metric cannot carry is refused") { + // A general triclinic metric has no monoclinic-C setting, and an F-centred cubic lattice has no + // hexagonal-P one (its hexagonal description is R-centred). + const CrystalLattice tri(41.0, 47.0, 53.0, 71.0, 83.0, 97.0); + CHECK_FALSE(LatticeSearchForClass(tri, gemmi::CrystalSystem::Monoclinic, 'C').has_value()); + + const double a = 60.0; + const auto cubic_f = CrystalLattice(a, a, a, 90, 90, 90).ToPrimitive('F'); + CHECK(LatticeSearch(cubic_f, 1e-6).centering == 'F'); + CHECK_FALSE(LatticeSearchForClass(cubic_f, gemmi::CrystalSystem::Hexagonal, 'P').has_value()); + // ... but its rhombohedral setting is there, which is what makes the refusal above a real answer + // rather than an artefact of the filter. + const auto rhomb = LatticeSearchForClass(cubic_f, gemmi::CrystalSystem::Trigonal, 'R'); + REQUIRE(rhomb.has_value()); + CHECK(rhomb->centering == 'R'); +} + +TEST_CASE("LatticeSearchForClass - asking for what the plain search found returns the same setting") { + const double a = 40.0; + const auto L = CrystalLattice(a, a, a, 90, 90, 90).ToPrimitive('I'); + const auto plain = LatticeSearch(L, 1e-6); + const auto filtered = LatticeSearchForClass(L, plain.system, plain.centering, 1e-6); + REQUIRE(filtered.has_value()); + CHECK(filtered->niggli_class == plain.niggli_class); + check_uc(filtered->conventional.GetUnitCell(), a, a, a, 90, 90, 90, 1e-4, 1e-4); +}