From 1b6831a2d3770e6d0708160b4c2eca9c2e80731c Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Fri, 10 Jul 2026 15:20:10 +0200 Subject: [PATCH] rugnux: gate the P1 space-group-search scaling at >= 1 The de-novo space-group determination (point-group intensity correlation plus the second-moment / twinning test) ran on the full resolution range, including noise-dominated high-res shells. There a single artifact - a zinger the box sum integrates whole, or a weak-reflection profile-fit runaway - has an astronomical resolution-normalised E (I / _shell with _shell ~ 0) and wrecks the normalised statistics, so the search collapses to P1. Re-scale the P1 search pass with a >= 1 cut over thin resolution shells (POINTLESS-style: determine symmetry where there is signal). Only this pass is cut; the final in-symmetry merge keeps the full range, so no output resolution is lost. A manual --scaling-high-resolution, if coarser, wins. Battery of 24 rotation crystals: box-sum lysoC_14 P1 -> P41212; and with the default gaussian integrator three previously mis-assigned crystals now get the correct space group (cytC_10 155->152, EP_cs_01-24 5->4, Ins_I_3 196->197), with no space-group regressions. Co-Authored-By: Claude Opus 4.8 (1M context) --- .../scale_merge/RotationScaleMerge.h | 4 ++ rugnux/Rugnux.cpp | 38 ++++++++++++++++++- 2 files changed, 41 insertions(+), 1 deletion(-) diff --git a/image_analysis/scale_merge/RotationScaleMerge.h b/image_analysis/scale_merge/RotationScaleMerge.h index f40b046c..9e4294d3 100644 --- a/image_analysis/scale_merge/RotationScaleMerge.h +++ b/image_analysis/scale_merge/RotationScaleMerge.h @@ -65,6 +65,10 @@ public: // masked_ice_rings: rings (indices into ICE_RING_RES_A) to drop from the final merge; empty = none. Result Run(bool for_search, const std::vector &masked_ice_rings = {}); + // Override the high-resolution cut for the next Run() - used to gate the de-novo P1 search pass at + // >= 1 without cutting the final in-symmetry merge. Reset to the manual limit afterwards. + void SetDMinLimit(std::optional d_min_A) { d_min_limit = d_min_A; } + private: // One integrated observation - a per-frame partial during scaling/combine, or a combined full during // scale-fulls/merge. Flat (not nested per image); a POD so the arrays translate straight to CUDA. diff --git a/rugnux/Rugnux.cpp b/rugnux/Rugnux.cpp index 3ecdacd3..1bc5004b 100644 --- a/rugnux/Rugnux.cpp +++ b/rugnux/Rugnux.cpp @@ -538,11 +538,47 @@ ProcessResult Rugnux::Run(RugnuxObserver *observer) { const auto initial_sg = experiment_.GetGemmiSpaceGroup(); auto sm = scale_and_merge(initial_sg ? initial_sg->short_name() : "P1", !initial_sg.has_value()); + // For a de-novo search, gate the P1 scaling+merge (the pass that feeds the space-group search) at + // >= 1 over thin resolution shells: noise-dominated high-res shells otherwise mislead + // the determination, where a single artifact (a zinger the box sum integrates whole, or a + // weak-reflection profile-fit runaway) has an astronomical resolution-normalised E and wrecks the + // intensity-correlation / second-moment statistics. The final in-symmetry merge keeps the full + // range (only this P1 search pass is cut); a manual --scaling-high-resolution, if coarser, wins. + double d_min_search = 0.0; + if (!initial_sg.has_value()) { + std::vector> rs; // (d, I/sigma) over the full P1 merge + rs.reserve(sm.merged.size()); + for (const auto &m : sm.merged) + if (std::isfinite(m.I) && std::isfinite(m.sigma) && m.sigma > 0.0f && + std::isfinite(m.d) && m.d > 0.0f) + rs.emplace_back(m.d, m.I / m.sigma); + if (rs.size() >= 400) { + std::sort(rs.begin(), rs.end(), + [](const auto &a, const auto &b) { return a.first > b.first; }); // low -> high res + const int bins = std::clamp(static_cast(rs.size() / 100), 1, 40); + const size_t per = (rs.size() + bins - 1) / static_cast(bins); + for (size_t b = 0; b * per < rs.size(); ++b) { + const size_t lo = b * per, hi = std::min(rs.size(), lo + per); + double sum = 0.0; + for (size_t j = lo; j < hi; ++j) sum += rs[j].second; + if (sum / static_cast(hi - lo) < 1.0) { d_min_search = rs[lo].first; break; } + } + } + const double d_manual = experiment_.GetScalingSettings().GetHighResolutionLimit_A().value_or(0.0); + if (rsm && d_min_search > d_manual) { + logger.Info("Space-group search: re-scaling the P1 pass at {:.2f} A ( >= 1)", d_min_search); + rsm->SetDMinLimit(d_min_search); + sm = scale_and_merge("P1, space-group search", true); + rsm->SetDMinLimit(d_manual > 0.0 ? std::optional(d_manual) : std::nullopt); + } + } + std::ostringstream stats_text; if (!experiment_.GetGemmiSpaceGroup().has_value()) { SearchSpaceGroupOptions sg_opts; sg_opts.merge_friedel = experiment_.GetScalingSettings().GetMergeFriedel(); - sg_opts.d_min_limit_A = experiment_.GetScalingSettings().GetHighResolutionLimit_A().value_or(0.0); + sg_opts.d_min_limit_A = std::max( + d_min_search, experiment_.GetScalingSettings().GetHighResolutionLimit_A().value_or(0.0)); // Constrain the search to subgroups of the lattice (metric) symmetry found by rotation // indexing. Centering is not constrained here - it is determined from the absences. if (end_msg.rotation_lattice_type.has_value())