From 093ac1612c84dc5aaafd5a3541ce2ef0885d4b96 Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Mon, 26 Oct 2026 20:20:08 +0100 Subject: [PATCH] Space-group search and frame disposition: three host-side waits shortened The first arm's space-group search now runs beside the wait for the all-observation merge and the search started on it, instead of after both. Classifying a reflection for a candidate group asks gemmi twice which operators fix it (once for the absence, once for the epsilon factor); one pass over the operators answers both, with the absence read off a fixing operator's own translation, which is what gemmi's test reduces to once the centring has been checked. The delta-CC1/2 frame disposition builds its per-reflection totals on all threads, each taking a range of reflections and walking every observation so the sums keep their order, and measures up to sixteen candidate chains side by side with a worker's scratch made on first use. Same numbers, same bits; only the waiting changed. Co-Authored-By: Claude Fable 5.1 Claude-Session: https://claude.ai/code/session_015CeEGgpKdfkCTJb3YNMJcu --- .../scale_merge/RotationScaleMerge.cpp | 39 ++++++++++++------- .../scale_merge/SearchSpaceGroup.cpp | 16 +++++++- rugnux/RugnuxScaleMerge.cpp | 8 +++- 3 files changed, 44 insertions(+), 19 deletions(-) diff --git a/image_analysis/scale_merge/RotationScaleMerge.cpp b/image_analysis/scale_merge/RotationScaleMerge.cpp index affc3f593..b1b9dd74d 100644 --- a/image_analysis/scale_merge/RotationScaleMerge.cpp +++ b/image_analysis/scale_merge/RotationScaleMerge.cpp @@ -3214,19 +3214,27 @@ void RotationScaleMerge::MeasureBatchDeltaCCHalf(int n_groups) { std::vector d; std::vector touched; }; - const int n_workers = static_cast(std::clamp(nthreads, 1, 8)); + // A worker's scratch is as long as the merge, so it is made on first use, by the worker: the + // workers a measurement does not need cost nothing, and the ones it does zero theirs side by side. + const int n_workers = static_cast(std::clamp(std::min(nthreads, n_batch), 1, 16)); std::vector scratch(n_workers); - for (auto &sc : scratch) - sc.d.assign(n_groups, Sums{}); + auto scratch_of = [&](int t) -> RangeScratch & { + if (scratch[t].d.empty()) scratch[t].d.assign(n_groups, Sums{}); + return scratch[t]; + }; + // Each worker takes a range of the reflections and walks every observation for them, so a + // reflection's sums are added in observation order whatever the split: the one serial sum. auto build_totals = [&]() { - std::fill(tot.begin(), tot.end(), Sums{}); - for (int i = 0; i < n_full; ++i) { - const int g = obs_group[i]; - if (g < 0 || frame_rejected[obs_frame[i]]) continue; - Sums &t = tot[g]; - t.n += 1; t.s1 += obs_I[i]; t.s2 += obs_I[i] * obs_I[i]; - t.w += obs_info[i]; - } + ParallelChunks(n_groups, nthreads, [&](int g0, int g1) { + std::fill(tot.begin() + g0, tot.begin() + g1, Sums{}); + for (int i = 0; i < n_full; ++i) { + const int g = obs_group[i]; + if (g < g0 || g >= g1 || frame_rejected[obs_frame[i]]) continue; + Sums &t = tot[g]; + t.n += 1; t.s1 += obs_I[i]; t.s2 += obs_I[i] * obs_I[i]; + t.w += obs_info[i]; + } + }); }; // delta-CC1/2 of the frames [f0, f1]. delta_cc is the difference in CC units, which is what a user @@ -3289,9 +3297,10 @@ void RotationScaleMerge::MeasureBatchDeltaCCHalf(int n_groups) { return r; }; auto measure = [&](int f0, int f1) { - accumulate(scratch[0], f0, f1); - const RangeDelta r = evaluate(scratch[0]); - clear(scratch[0]); + RangeScratch &sc = scratch_of(0); + accumulate(sc, f0, f1); + const RangeDelta r = evaluate(sc); + clear(sc); return r; }; // Every range in `ranges`, side by side. Ranges that start on the same frame are one chain, measured @@ -3311,7 +3320,7 @@ void RotationScaleMerge::MeasureBatchDeltaCCHalf(int n_groups) { std::vector out(ranges.size()); std::atomic next_chain{0}; parallel_detail::RunTasks(std::min(n_workers, n_chain), [&](int t) { - RangeScratch &sc = scratch[t]; + RangeScratch &sc = scratch_of(t); for (int c = next_chain++; c < n_chain; c = next_chain++) { int added_to = ranges[sorted[chain_start[c]]].first - 1; for (size_t k = chain_start[c]; k < chain_start[c + 1]; ++k) { diff --git a/image_analysis/scale_merge/SearchSpaceGroup.cpp b/image_analysis/scale_merge/SearchSpaceGroup.cpp index 6cf37304d..92b7a3bb1 100644 --- a/image_analysis/scale_merge/SearchSpaceGroup.cpp +++ b/image_analysis/scale_merge/SearchSpaceGroup.cpp @@ -144,9 +144,21 @@ namespace { if (CenteringAbsent(gops, hkl)) return REFL_CENTRING_ABSENT; const uint8_t zone = static_cast((MirrorZoneOf(mirrors, hkl) + 1) << 4); - if (gops.is_systematically_absent(hkl)) + // gemmi's is_systematically_absent and epsilon_factor_without_centering in one pass over the + // operators: both ask which of them fix hkl. The absence test adds each centring vector to a + // fixing operator's translation, but no centring vector shifts the phase of this reflection + // (above), so the operator's own translation decides, as it does for one with none. + const gemmi::Op::Miller denh = {{gemmi::Op::DEN * hkl[0], gemmi::Op::DEN * hkl[1], gemmi::Op::DEN * hkl[2]}}; + int epsilon = 0; + bool absent = false; + for (const gemmi::Op& op : gops.sym_ops) + if (op.apply_to_hkl_without_division(hkl) == denh) { + ++epsilon; + absent |= gemmi::GroupOps::has_phase_shift(op.tran, hkl); + } + if (absent) return REFL_ABSENT | zone; - return REFL_PRESENT | zone | (gops.epsilon_factor_without_centering(hkl) > 1 ? REFL_ON_AXIS : 0); + return REFL_PRESENT | zone | (epsilon > 1 ? REFL_ON_AXIS : 0); } // The zone a mirror fixes, named the way a crystallographer names it - "h0l" for the mirror diff --git a/rugnux/RugnuxScaleMerge.cpp b/rugnux/RugnuxScaleMerge.cpp index 48b94dd89..51d157475 100644 --- a/rugnux/RugnuxScaleMerge.cpp +++ b/rugnux/RugnuxScaleMerge.cpp @@ -1916,6 +1916,11 @@ bool Rugnux::ScaleMergeAndSymmetry(PipelineLocals &p) { // else. Whether the all-observation arm is wanted at all is known only from the first arm // (below), so nothing of it reaches the run before then. A merge that failed beside the // others is made again where it is wanted, in series, as it would have been without it. + // The first arm's search runs beside the wait for that merge and the search on it. + TimingMark first_arm_mark(logger, "SearchSpaceGroup (first arm)"); + auto first_search_ahead = std::async(std::launch::async, [&merged = sm.merged, sg_opts] { + return SearchSpaceGroup(merged, sg_opts); + }); std::optional all_merge; std::future all_search_ahead; if (rsm_all && !all_arm_ahead.valid()) @@ -1937,8 +1942,7 @@ bool Rugnux::ScaleMergeAndSymmetry(PipelineLocals &p) { return SearchSpaceGroup(merged, all_opts); }); } - TimingMark first_arm_mark(logger, "SearchSpaceGroup (first arm)"); - auto sg_search = SearchSpaceGroup(sm.merged, sg_opts); + auto sg_search = first_search_ahead.get(); first_arm_mark.Stop(); // The Lorentz-filtered arm's merge, kept past the point where the all-observation one