From 8c288c8697b9caa216a540feb433a439266f6a49 Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Thu, 9 Jul 2026 21:14:18 +0200 Subject: [PATCH] rugnux: recover genuine high symmetry in de-novo space-group search The POINTLESS-style symmetry search discarded a real point group on two hard rotation-test crystals: - Ins_I_3 (I23): one of the eight cubic 3-fold operators scored CC 0.483, just under the 0.5 per-operator floor, so I23 fell to I222 and merging collapsed (3-fold orientation ambiguity; CC1/2 1%, R-meas 558%). Gate A now scores each CONJUGACY CLASS of a point group's rotations by its mean CC - conjugate operators are symmetry-equivalent and must stand or fall together - so the cubic 3-fold class (mean 0.57) is confirmed. -> I23, CC1/2 99.9%. Degenerate (few-pair / NaN) class members are skipped rather than vetoing the class. - Benas_3/7 (F432): the chi^2 self-consistency gate rejected the true cubic group because the weak-data error model is badly miscalibrated (reduced chi^2 ~150), which inflates the chi^2 ratio with point-group order for genuine high symmetry too. The ratio bound is now widened by log10(chi2_ref) so a broken error model does not veto real symmetry, while a false operator whose ratio is far worse than its subgroup's is still caught. The bound stays a per-candidate chi^2 test - chi^2 still arbitrates and is never bypassed. Verified: F432 (chi2_ref 151, ratio 3.2 < bound 4.2) accepted; a calibrated pseudo/twin R32 (Ins_H_3: chi2_ref 2.6, ratio 4.2 > bound 2.4) still rejected, stays R3. SearchSpaceGroup.{cpp,h} only. Validated on the 24-crystal rotation-test battery: exactly these two space groups change (both now match the reference); zero regressions - every pseudo-symmetry P1 case and every correct assignment unchanged. Co-Authored-By: Claude Opus 4.8 (1M context) --- .../scale_merge/SearchSpaceGroup.cpp | 81 +++++++++++++++---- image_analysis/scale_merge/SearchSpaceGroup.h | 3 + 2 files changed, 69 insertions(+), 15 deletions(-) diff --git a/image_analysis/scale_merge/SearchSpaceGroup.cpp b/image_analysis/scale_merge/SearchSpaceGroup.cpp index ec434941..e407cd56 100644 --- a/image_analysis/scale_merge/SearchSpaceGroup.cpp +++ b/image_analysis/scale_merge/SearchSpaceGroup.cpp @@ -271,6 +271,54 @@ SearchSpaceGroupResult SearchSpaceGroup( return op_cache.emplace(rk, score_operator(op)).first->second; }; + // Conjugate rotations (symmetry-equivalent within the point group) relate symmetry-equivalent + // reflection sets, so on real data their CCs cluster; a noisy crystal can push one class member + // below min_operator_cc while the class is unmistakably present (e.g. one cubic 3-fold at 0.48 + // among siblings at 0.53-0.66). Judge each conjugacy class by its mean CC, not its weakest + // member, so a genuine high-symmetry point group is not lost to one marginal operator. chi2_under + // (below) remains the safety net against a truly false promotion. Returns {all classes present, + // worst class-mean CC}. + auto point_group_present = [&](const std::vector& rots) -> std::pair { + const size_t m = rots.size(); + std::vector cls(m, -1); + int n_cls = 0; + for (size_t i = 0; i < m; ++i) { + if (cls[i] >= 0) + continue; + cls[i] = n_cls; + for (size_t j = i + 1; j < m; ++j) + if (cls[j] < 0) + for (const auto& p : rots) + if ((p * rots[i] * p.inverse()).rot == rots[j].rot) { + cls[j] = n_cls; + break; + } + ++n_cls; + } + bool ok = true; + double worst_mean = 1.0; + for (int c = 0; c < n_cls; ++c) { + // Average over the class members that actually have enough pairs to score; a single + // low-multiplicity / degenerate (NaN) operator in an otherwise strong class is skipped, + // not allowed to veto the class. The class must still have at least one scored member. + double sum_cc = 0.0; + int n_valid = 0; + for (size_t i = 0; i < m; ++i) + if (cls[i] == c) { + const auto& s = operator_score(rots[i]); + if (s.n_pairs >= opt.min_pairs_per_operator && std::isfinite(s.cc)) { + sum_cc += s.cc; + ++n_valid; + } + } + const double mean_cc = n_valid > 0 ? sum_cc / n_valid : 0.0; + worst_mean = std::min(worst_mean, mean_cc); + if (n_valid == 0 || mean_cc < opt.min_operator_cc) + ok = false; + } + return {ok, worst_mean}; + }; + std::optional holohedry; if (opt.lattice_system.has_value()) holohedry = HolohedryRotationSet(opt.lattice_system.value()); @@ -319,41 +367,44 @@ SearchSpaceGroupResult SearchSpaceGroup( }; // Operator-CC-confirmed candidates, each with its merge chi^2; chi2_ref = the most consistent. - struct PGCand { const PointGroupInfo* pg; int order; double min_cc; double chi2; }; + struct PGCand { const PointGroupInfo* pg; int order; double min_class_cc; double chi2; }; std::vector pg_cands; double chi2_ref = std::numeric_limits::infinity(); for (const auto& pg : point_groups) { - bool all_present = true; - double min_cc = pg.rotations.empty() ? 1.0 : std::numeric_limits::infinity(); - for (const auto& op : pg.rotations) { - const auto& s = operator_score(op); - all_present = all_present && s.present; - min_cc = std::min(min_cc, s.cc); - } - if (!all_present) + const auto [present, min_class_cc] = point_group_present(pg.rotations); + if (!present) continue; const double ch = pg.rotations.empty() ? std::numeric_limits::quiet_NaN() : chi2_under(pg.rotations); - pg_cands.push_back({&pg, static_cast(pg.rotations.size()) + 1, min_cc, ch}); + pg_cands.push_back({&pg, static_cast(pg.rotations.size()) + 1, min_class_cc, ch}); if (!pg.rotations.empty() && std::isfinite(ch)) chi2_ref = std::min(chi2_ref, ch); } // Choose the largest point group that is both operator-confirmed AND self-consistent (its merge - // chi^2 is not inflated past max_merge_chi2_ratio x the most-consistent candidate; ties -> higher - // min CC). Identity (no operators) is always consistent, so it stays the P1 fallback. + // chi^2 is not inflated past the miscalibration-widened bound below; ties -> higher min class CC). + // Identity (no operators) is always consistent, so it stays the P1 fallback. const PointGroupInfo* best_pg = nullptr; int best_pg_order = 0; double best_pg_min_cc = -2.0; for (const auto& c : pg_cands) { + // Widen the chi^2 ratio bound by the error-model miscalibration. When even the best candidate's + // reduced chi^2 (chi2_ref) is far above 1 (weak data / uncorrected decay), the ratio grows with + // point-group order for genuine high symmetry too, so the fixed bound wrongly rejects it. Adding + // log10(chi2_ref) lets the bound exceed max_merge_chi2_ratio by ~the order of magnitude of the + // miscalibration - enough to keep the true group (F432 weak data: chi2_ref~150, ratio 3.2 vs + // bound ~4.2) while STILL rejecting a false operator whose ratio is far worse than its subgroup's + // (a twin/pseudo R32 on calibrated data: chi2_ref~2.6, ratio 4.2 vs bound ~2.4). The bound stays + // a per-candidate chi^2 test, so chi^2 still arbitrates between candidates - it is never bypassed. + const double allowed_ratio = opt.max_merge_chi2_ratio + std::log10(std::max(1.0, chi2_ref)); const bool consistent = c.pg->rotations.empty() || !std::isfinite(c.chi2) || - !std::isfinite(chi2_ref) || c.chi2 <= chi2_ref * opt.max_merge_chi2_ratio; + !std::isfinite(chi2_ref) || c.chi2 <= chi2_ref * allowed_ratio; if (!consistent) continue; - if (c.order > best_pg_order || (c.order == best_pg_order && c.min_cc > best_pg_min_cc)) { + if (c.order > best_pg_order || (c.order == best_pg_order && c.min_class_cc > best_pg_min_cc)) { best_pg = c.pg; best_pg_order = c.order; - best_pg_min_cc = c.min_cc; + best_pg_min_cc = c.min_class_cc; } } diff --git a/image_analysis/scale_merge/SearchSpaceGroup.h b/image_analysis/scale_merge/SearchSpaceGroup.h index c512525f..021e38ea 100644 --- a/image_analysis/scale_merge/SearchSpaceGroup.h +++ b/image_analysis/scale_merge/SearchSpaceGroup.h @@ -79,6 +79,9 @@ struct SearchSpaceGroupOptions { // sigma and chi^2 blows up; a real one leaves it ~flat even when the operator CC is only moderate. // 2.0 separates the test set: a true point group's chi^2 stays within ~1.6x the best subgroup // (the cubic merge of a noisy crystal is the worst real case), while a false operator gives >=2.3x. + // This bound is the calibrated-error baseline; SearchSpaceGroup widens it by log10(chi2_ref) so a + // badly miscalibrated error model (weak data / uncorrected decay, where the ratio grows with + // point-group order for real symmetry too) does not spuriously reject a genuine high-symmetry group. double max_merge_chi2_ratio = 2.0; // --- Stage B: space group (screw axes / centering) ---