Merge branch 'perf-scale-c1' into rc174

docs/CPU_DATA_ANALYSIS_INTEGRATION.md

@@ -384,7 +384,7 @@ The combine groups each reflection's partials into rocking events (contiguous ru
 
 The fulls are then re-scaled in the XDS sense — a per-image scale refit directly on the complete reflections under the unity partiality model — and merged (§10.4).
 
-Before the combine a per-frame scale can also be fitted on the partials themselves. That needs a frame to hold many rocking events caught at different points of their curves: within one rocking curve a change of scale and an error of the partiality model are the same thing, and on a finely sliced sparse sweep the fit takes one for the other. The rocking events per frame are the prior (the partials are scaled from 50 events per frame), but the counts of small-molecule sweeps (3–43) and proteins (5–900) overlap. So the first merge of a run that is not a space-group search is made **both ways**, with the partials scaled and with the scale taken from the fulls alone, and the prior stands unless its merge has no resolved error model ($b$ not resolved from zero: the strong equivalents do not agree to within a measurable systematic error) while the other merge has one. The two ISa values are deliberately not compared beyond that: the partiality-model error a partial scale takes up is shared by symmetry mates measured at the same rocking geometry, so their agreement cannot see it — a small-molecule sweep with six events per frame read ISa 10.5 with its partials scaled against 8.7 from the fulls alone, and refined to $R_1$ 0.105 against 0.062. The log states the choice and both arms' ISa. Because every merged observation is now a counting-statistics-limited full rather than a partiality-divided slice, the error model reaches a far higher asymptotic $I/\sigma$.
+Before the combine a per-frame scale can also be fitted on the partials themselves. That needs a frame to hold many rocking events caught at different points of their curves: within one rocking curve a change of scale and an error of the partiality model are the same thing, and on a finely sliced sparse sweep the fit takes one for the other. The rocking events per frame are the prior (the partials are scaled from 50 events per frame), but the counts of small-molecule sweeps (3–43) and proteins (5–900) overlap. So the first merge of a run that is not a space-group search can be made **both ways**, with the partials scaled and with the scale taken from the fulls alone, and the prior stands unless its merge has no resolved error model ($b$ not resolved from zero: the strong equivalents do not agree to within a measurable systematic error) while the other merge has one. Where the prior is to scale the partials and that merge resolves its error model, the other cannot change the choice and is not made. The two ISa values are deliberately not compared beyond that: the partiality-model error a partial scale takes up is shared by symmetry mates measured at the same rocking geometry, so their agreement cannot see it — a small-molecule sweep with six events per frame read ISa 10.5 with its partials scaled against 8.7 from the fulls alone, and refined to $R_1$ 0.105 against 0.062. The log states the choice and the ISa of every arm that was made. Because every merged observation is now a counting-statistics-limited full rather than a partiality-divided slice, the error model reaches a far higher asymptotic $I/\sigma$.
 
 How smooth that scale is over the rotation is left to the data rather than to a fixed window. Each round fits every frame on its own fulls against the current reference, giving a scale $G_f$ and its information $D_f=\sum w^2c^2$; the scale is then the curve $x=\log G$ minimising $\sum_f J_f\,(x_f-y_f)^2+\lambda\sum_f(\Delta^2 x)_f^2$, with $y_f$ the frame's own fit in log scale and $J_f$ its information carried there — a penalised (Whittaker–Eilers) smoother, solved as a five-band linear system. $\lambda$ is chosen by cross-validation: blocks of frames one rocking curve wide are left out in turn and predicted from the curve through the rest (neighbours closer than a rocking curve share their measurement, since a full sums those frames). A frame of hundreds of fulls is then followed frame by frame, a frame of two or three is carried by its neighbours, a stretch with none is bridged by a straight line, and a scale that falls a hundredfold over a few degrees — an absorbing crystal turning edge-on — is followed where a window would average across it. Once the curve settles, one free fit is shrunk toward it frame by frame by how much of each frame's deviation its neighbour shares (the lag-1 covariance), which hands back a real per-frame systematic and discards fit noise.
 

image_analysis/scale_merge/RotationScaleMerge.cpp

@@ -233,7 +233,7 @@ RotationScaleMerge::RotationScaleMerge(const DiffractionExperiment &experiment,
                                        size_t nthreads, Logger &logger,
                                        std::string observation_dump_path)
     : x(experiment), partials_out(partial_outcomes), reference_cell(std::move(reference_cell)),
-      nthreads(nthreads == 0 ? std::thread::hardware_concurrency() : nthreads), logger(logger),
+      nthreads(nthreads == 0 ? std::thread::hardware_concurrency() : nthreads), logger(&logger),
       observation_dump_path(std::move(observation_dump_path)) {
     const auto s = x.GetScalingSettings();
     min_partiality = s.GetMinPartiality();
@@ -509,7 +509,7 @@ void RotationScaleMerge::Ingest() {
     }
     BuildInRangeObservations(keys, static_cast<int>(total));
     rawrun_group.assign(rawrun_start.size(), -1);
-    logger.Info("RotationScaleMerge: ingested {} partial observations from {} frames ({} distinct hkl)",
+    logger->Info("RotationScaleMerge: ingested {} partial observations from {} frames ({} distinct hkl)",
                 n_partials_obs, n_frames, rawrun_start.size());
 
     // The partiality the predictor gave each partial, which its corr was formed with (image_scale_corr
@@ -616,7 +616,7 @@ void RotationScaleMerge::Ingest() {
                          rawrun_start.data(), rawrun_count.data(),
                          rawrun_h.data(), rawrun_k.data(), rawrun_l.data());
         gpu_->SetFrameCellOk(frame_cell_ok.data());
-        logger.Info("RotationScaleMerge: GPU scaling + combine + scale-fulls + merge active");
+        logger->Info("RotationScaleMerge: GPU scaling + combine + scale-fulls + merge active");
     }
     // The device now holds every per-obs field the resident pipeline reads, so the ingest arrays'
     // one remaining reader is RockingEventFrames - which on this path sees only ingest-time values
@@ -875,7 +875,7 @@ void RotationScaleMerge::BuildInRangeObservations(std::unique_ptr<SortKey[]> &ke
     ingest_keep = std::move(keep);
 
     if (n_keep < n_obs)
-        logger.Info("RotationScaleMerge: dropped {} of {} observations ({} of {} distinct hkl) outside the "
+        logger->Info("RotationScaleMerge: dropped {} of {} observations ({} of {} distinct hkl) outside the "
                     "requested resolution range - nothing downstream could have merged them",
                     n_obs - n_keep, n_obs, n_run - n_keep_run, n_run);
 }
@@ -1065,7 +1065,7 @@ void RotationScaleMerge::SmoothGeometry() {
         changed += local;
     });
 
-    logger.Info("Smoothed per-frame geometry over +-{} frames (chosen by cross-validation); "
+    logger->Info("Smoothed per-frame geometry over +-{} frames (chosen by cross-validation); "
                 "recomputed delta_phi for {} of {} partials", best_half, changed.load(),
                 n_partials_obs);
 }
@@ -1127,7 +1127,7 @@ void RotationScaleMerge::SmoothMosaicityAndPartiality() {
             n_events += l_ev;
             n_partials += l_pa;
         });
-        logger.Info("One exact-Bragg angle per rocking event: {} events, {} partials relaid",
+        logger->Info("One exact-Bragg angle per rocking event: {} events, {} partials relaid",
                     n_events.load(), n_partials.load());
     }
 
@@ -1202,7 +1202,7 @@ void RotationScaleMerge::SmoothMosaicityAndPartiality() {
                 ingest_partiality[i] = RotationPartiality(ingest_delta_phi[i], r.zeta, mos, wedge);
         });
     });
-    logger.Info("Recomputed partiality from frame-order-smoothed mosaicity");
+    logger->Info("Recomputed partiality from frame-order-smoothed mosaicity");
 }
 
 int RotationScaleMerge::ComputeAsuGroups(const HKLKeyGenerator &keygen) {
@@ -1230,14 +1230,17 @@ int RotationScaleMerge::ComputeAsuGroups(const HKLKeyGenerator &keygen) {
     // The key travels with the run it belongs to, so the sort and the run-detection below read the
     // key they are standing on instead of chasing it back through `key` - an indirect compare is a
     // cache miss per comparison, and there are a few tens of millions of them. Which of two runs
-    // with the SAME key ends up first is left to the sort, and cannot matter: the packed key holds
-    // h, k, l and the hand exactly, so every run in a tie reduces to the same ASU reflection.
+    // with the SAME key ends up first cannot matter: the packed key holds h, k, l and the hand
+    // exactly, so every run in a tie reduces to the same ASU reflection. The tie is still broken on
+    // the run, so the order is a total one and the sort can go on all threads (ParallelSort).
     struct RunKey { uint64_t key; int32_t run; };
     std::vector<RunKey> sorted;
     sorted.reserve(n_run);
     for (int r = 0; r < n_run; ++r)
         if (eligible[r]) sorted.push_back({key[r], r});
-    std::sort(sorted.begin(), sorted.end(), [](const RunKey &a, const RunKey &b) { return a.key < b.key; });
+    ParallelSort(sorted.begin(), sorted.end(), nthreads, [](const RunKey &a, const RunKey &b) {
+        return a.key < b.key || (a.key == b.key && a.run < b.run);
+    });
 
     // Where a group starts is a property of the sorted keys alone, so mark the starts on all threads and
     // number them with a prefix - the same ids the append-as-you-go walk handed out. The gemmi ASU
@@ -1624,7 +1627,7 @@ void RotationScaleMerge::ShrinkToRestrained(const char *what, const std::vector<
         kept_sum += kept;
         ++kept_n;
     }
-    logger.Info("Per-frame scaling ({}): the scales carry {:.0f}% of what each frame's own "
+    logger->Info("Per-frame scaling ({}): the scales carry {:.0f}% of what each frame's own "
                 "observations say beyond its neighbourhood's - that deviation is {:.2f}% rms, of "
                 "which {:.2f}% is shared with the next frame and the rest is the fit's own noise",
                 what, 100.0 * (kept_n ? kept_sum / kept_n : 0.0), 100.0 * std::sqrt(var_d),
@@ -1721,7 +1724,7 @@ RotationScaleMerge::ScalingLoopOutcome RotationScaleMerge::RunScalingLoop(
                     ++n;
                 }
             out.step = sw > 0.0 ? std::sqrt(ss / sw) : 0.0;
-            logger.Debug("Per-frame scaling ({}): iteration {}, rms |dlogG| {:.2e} over {} frames, "
+            logger->Debug("Per-frame scaling ({}): iteration {}, rms |dlogG| {:.2e} over {} frames, "
                          "gauge {:.4f}", what, it + 1, out.step, n, g_ref);
             if (out.step < SCALING_TOLERANCE) {
                 out.converged = true;
@@ -1745,15 +1748,15 @@ RotationScaleMerge::ScalingLoopOutcome RotationScaleMerge::RunScalingLoop(
 
 void RotationScaleMerge::LogScalingOutcome(const char *what, const ScalingLoopOutcome &out) const {
     if (out.converged)
-        logger.Info("Per-frame scaling ({}): settled after {} iterations (rms |dlogG| {:.1e})", what,
+        logger->Info("Per-frame scaling ({}): settled after {} iterations (rms |dlogG| {:.1e})", what,
                     out.iterations, out.step);
     else if (out.stalled)
-        logger.Warning("Per-frame scaling ({}): stopped after {} iterations - the scales stopped settling "
+        logger->Warning("Per-frame scaling ({}): stopped after {} iterations - the scales stopped settling "
                        "(rms |dlogG| {:.1e}, not below {:.1e} for five rounds); more rounds would only walk "
                        "them, so everything read off this merge is read off an unsettled state", what,
                        out.iterations, out.step, out.best_step);
     else
-        logger.Warning("Per-frame scaling ({}): NOT settled after the cap of {} iterations - the scales "
+        logger->Warning("Per-frame scaling ({}): NOT settled after the cap of {} iterations - the scales "
                        "were still moving by rms |dlogG| {:.1e} (tolerance {:.0e}); everything read off "
                        "this merge is read off an unsettled state", what, out.iterations, out.step,
                        SCALING_TOLERANCE);
@@ -1823,7 +1826,7 @@ RotationScaleMerge::ScalingLoopOutcome RotationScaleMerge::RunFullsScalingLoop(
             }
         g = g_new;
         out.step = sw > 0.0 ? std::sqrt(ss / sw) : 0.0;
-        logger.Debug("Per-frame scaling (fulls): iteration {}, rms |dlogG| {:.2e}, lambda {:.3g}, gauge {:.4f}",
+        logger->Debug("Per-frame scaling (fulls): iteration {}, rms |dlogG| {:.2e}, lambda {:.3g}, gauge {:.4f}",
                      it + 1, out.step, lambda, g_ref);
         if (out.step < SCALING_TOLERANCE) {
             out.converged = true;
@@ -1863,7 +1866,7 @@ RotationScaleMerge::ScalingLoopOutcome RotationScaleMerge::RunFullsScalingLoop(
     for (double j : J)
         if (j > 0.0) positive.push_back(j);
     if (!positive.empty() && lambda > 0.0)
-        logger.Info("Per-frame scaling (fulls): smoothness chosen by cross-validation over the frames - "
+        logger->Info("Per-frame scaling (fulls): smoothness chosen by cross-validation over the frames - "
                     "a frame's scale is carried over about {:.0f} frames", std::pow(lambda / median_of(positive), 0.25));
     LogScalingOutcome("fulls", out);
     return out;
@@ -2015,7 +2018,7 @@ void RotationScaleMerge::RefineDecay(int n_groups) {
     const double slope = fit_slope(-1);
     const double total_delta_B = std::fabs(0.5 * slope * n_frames);
     if (total_delta_B < DECAY_MIN_DELTA_B) {
-        logger.Info("Decay correction: negligible radiation damage (total dB = {:.2f} A^2 < {:.1f}, skipped)",
+        logger->Info("Decay correction: negligible radiation damage (total dB = {:.2f} A^2 < {:.1f}, skipped)",
                     total_delta_B, DECAY_MIN_DELTA_B);
         return;
     }
@@ -2024,7 +2027,7 @@ void RotationScaleMerge::RefineDecay(int n_groups) {
     const double base = subset_disagreement(1, 0.0) + subset_disagreement(0, 0.0);
     const double gain = base - (subset_disagreement(1, fit_slope(0)) + subset_disagreement(0, fit_slope(1)));
     if (!(gain > CV_MIN_RELATIVE_GAIN * base)) {
-        logger.Info("Decay correction: not cross-validated (dB = {:.2f} A^2, held-out gain {:.1f}%, skipped)",
+        logger->Info("Decay correction: not cross-validated (dB = {:.2f} A^2, held-out gain {:.1f}%, skipped)",
                     total_delta_B, 100.0 * gain / std::max(base, 1e-30));
         return;
     }
@@ -2033,7 +2036,7 @@ void RotationScaleMerge::RefineDecay(int n_groups) {
             if (usable(fulls[i]))
                 fulls[i].corr = static_cast<float>(fulls[i].corr * decay_factor(fulls[i], slope));
     });
-    logger.Info("Decay correction: total relative-B = {:.2f} A^2 over run (dB/dframe = {:.2e}, cross-validated)",
+    logger->Info("Decay correction: total relative-B = {:.2f} A^2 over run (dB/dframe = {:.2e}, cross-validated)",
                 0.5 * slope * n_frames, 0.5 * slope);
 }
 
@@ -3081,7 +3084,7 @@ void RotationScaleMerge::MeasureBatchDeltaCCHalf(int n_groups) {
         int n_new = 0;
         for (int f = f0; f <= f1; ++f) n_new += frame_rejected[f] ? 0 : 1;
         if (n_rejected + n_new > max_rejected) {
-            logger.Warning("delta-CC1/2 stopped after {} of {} frames: removing the next stretch would "
+            logger->Warning("delta-CC1/2 stopped after {} of {} frames: removing the next stretch would "
                            "take out more than {:.0f}% of the sweep, and past that the merge is no "
                            "longer a reference the rest of the run can be judged against",
                            n_rejected, n_frames, 100.0 * DELTA_CC_MAX_REJECTED);
@@ -3134,7 +3137,7 @@ void RotationScaleMerge::MeasureBatchDeltaCCHalf(int n_groups) {
     }
 
     if (n_rejected > 0)
-        logger.Info("delta-CC1/2: {} of {} frames ({:.1f} deg) removed from the merge over {} stretch(es) "
+        logger->Info("delta-CC1/2: {} of {} frames ({:.1f} deg) removed from the merge over {} stretch(es) "
                     "- keeping them lowered CC1/2 by more than {:.0f} standard errors",
                     n_rejected, n_frames, n_rejected * osc_deg, rejected_ranges.size(), DELTA_CC_SIGMA);
 }
@@ -3279,7 +3282,7 @@ void RotationScaleMerge::RefineRelativeB(int n_groups) {
     const double bmin = *std::min_element(b_all.begin(), b_all.end());
     const double bmax = *std::max_element(b_all.begin(), b_all.end());
     if (bmax - bmin < RELATIVE_B_MIN_SPREAD) {
-        logger.Info("Relative-B: negligible variation ({} batches, peak-to-peak {:.2f} A^2 < {:.1f}, skipped)",
+        logger->Info("Relative-B: negligible variation ({} batches, peak-to-peak {:.2f} A^2 < {:.1f}, skipped)",
                     n_batch, bmax - bmin, RELATIVE_B_MIN_SPREAD);
         return;
     }
@@ -3291,7 +3294,7 @@ void RotationScaleMerge::RefineRelativeB(int n_groups) {
     const std::vector<double> b_odd = FitRelativeBCurve(t, n_groups, n_batch, frames_per_batch, 1);
     const double gain = base - (subset_disagreement(1, b_even) + subset_disagreement(0, b_odd));
     if (!(gain > CV_MIN_RELATIVE_GAIN * base)) {
-        logger.Info("Relative-B: not cross-validated ({} batches, peak-to-peak {:.2f} A^2, held-out gain {:.1f}%, skipped)",
+        logger->Info("Relative-B: not cross-validated ({} batches, peak-to-peak {:.2f} A^2, held-out gain {:.1f}%, skipped)",
                     n_batch, bmax - bmin, 100.0 * gain / std::max(base, 1e-30));
         return;
     }
@@ -3300,7 +3303,7 @@ void RotationScaleMerge::RefineRelativeB(int n_groups) {
             if (usable(fulls[i]))
                 fulls[i].corr = static_cast<float>(fulls[i].corr * b_factor(fulls[i], b_all));
     });
-    logger.Info("Relative-B: {} batches of ~{:.1f} deg, peak-to-peak {:.2f} A^2, held-out gain {:.1f}% (cross-validated)",
+    logger->Info("Relative-B: {} batches of ~{:.1f} deg, peak-to-peak {:.2f} A^2, held-out gain {:.1f}% (cross-validated)",
                 n_batch, relative_b_deg, bmax - bmin, 100.0 * gain / std::max(base, 1e-30));
 }
 
@@ -3951,7 +3954,7 @@ void RotationScaleMerge::ApplyCellSurface(const std::vector<int32_t> &cell, int
     }
     if (n_cc > 0) gain /= n_cc;
     if (!(gain > 0.0)) {
-        logger.Info("{} correction: not cross-validated (held-out half-set CC {:+.4f} per shell on Fisher's z, "
+        logger->Info("{} correction: not cross-validated (held-out half-set CC {:+.4f} per shell on Fisher's z, "
                     "skipped)", name, gain);
         return;
     }
@@ -3961,7 +3964,7 @@ void RotationScaleMerge::ApplyCellSurface(const std::vector<int32_t> &cell, int
             if (cell[i] >= 0)
                 fulls[i].corr = static_cast<float>(fulls[i].corr * A[cell[i]]);
     });
-    logger.Info("{} correction: cross-validated, held-out half-set CC {:+.4f} per shell on Fisher's z; {} in {} round(s); "
+    logger->Info("{} correction: cross-validated, held-out half-set CC {:+.4f} per shell on Fisher's z; {} in {} round(s); "
                 "resolution gauge removed a {:.2f}x ramp; {} of {} cells on the [0.25, 4] clamp",
                 name, gain, settled ? "settled" : "NOT settled",
                 rounds_used, gauge_ramp, n_clamped, ncell);
@@ -4063,7 +4066,7 @@ bool RotationScaleMerge::DropCollapsedScales(const std::vector<uint8_t> &fitted_
         ++n_dropped;
     }
     if (n_dropped > 0)
-        logger.Warning("Dropped {} frame(s) whose scale came out more than {:.0f}x below the run median "
+        logger->Warning("Dropped {} frame(s) whose scale came out more than {:.0f}x below the run median "
                        "(worst {:.0f}x) - a scale that small says the frame holds no diffraction, and "
                        "using it would amplify the frame's intensities, and its sigmas with them, "
                        "by the same factor",
@@ -4132,7 +4135,7 @@ bool RotationScaleMerge::DropCollapsedFullScales(bool from_staging) {
             ++n_dropped;
         }
     if (n_dropped > 0)
-        logger.Warning("Dropped {} frame(s) / {} full(s) whose scale came out more than {:.0f}x below "
+        logger->Warning("Dropped {} frame(s) / {} full(s) whose scale came out more than {:.0f}x below "
                        "the run median (worst {:.0f}x) - a scale that small says the frame holds no "
                        "diffraction",
                        n_frames_dropped, n_dropped, 1.0 / MIN_CREDIBLE_SCALE_RATIO, worst);
@@ -4338,7 +4341,7 @@ void RotationScaleMerge::Combine() {
     }
 
     SortFullsByFrame();
-    logger.Info("3D combine: {} fulls from {} partials, {} rocking events dropped for a saturated pixel",
+    logger->Info("3D combine: {} fulls from {} partials, {} rocking events dropped for a saturated pixel",
                 fulls.size(), n_used, fulls_dropped_overloaded);
 }
 
@@ -4507,8 +4510,8 @@ RotationScaleMerge::Result RotationScaleMerge::MergeAndStats(int n_groups, bool
     std::vector<uint8_t> obs_hand(fulls.size(), 0);        // per-full, 0 = I(+), 1 = I(-)
     std::vector<uint8_t> group_has_hands(n_groups, 0);
     if (merge_friedel) {
-        const HKLKeyGenerator anom_keygen(/*merge_friedel=*/false, x.GetSpaceGroupOrP1());
-        const gemmi::GroupOps gops = x.GetSpaceGroupOrP1().operations();
+        const HKLKeyGenerator anom_keygen(/*merge_friedel=*/false, Group());
+        const gemmi::GroupOps gops = Group().operations();
         ParallelChunks(n_groups, ThreadsForWork(n_groups, nthreads), [&](int lo, int hi) {
             for (int g = lo; g < hi; ++g)
                 group_has_hands[g] = gops.is_reflection_centric(
@@ -4680,7 +4683,7 @@ RotationScaleMerge::Result RotationScaleMerge::MergeAndStats(int n_groups, bool
             // set n_groups + p is Friedel pair p's, filled only where a hand of p needs it.
             std::vector<uint8_t> pair_needed;
             if (!merge_friedel) {
-                const HKLKeyGenerator anom_keygen(/*merge_friedel=*/false, x.GetSpaceGroupOrP1());
+                const HKLKeyGenerator anom_keygen(/*merge_friedel=*/false, Group());
                 std::unordered_map<uint64_t, int32_t> pair_id;
                 pair_id.reserve(static_cast<size_t>(n_groups) + 1);
                 pair_of_group.resize(n_groups);
@@ -4772,7 +4775,9 @@ RotationScaleMerge::Result RotationScaleMerge::MergeAndStats(int n_groups, bool
                 lv.emplace_back(std::fabs(std::log(I_corr / reject_median[g])), em_mean[g] * em_mean[g] / var);
             }
             if (!lv.empty()) {
-                std::sort(lv.begin(), lv.end());
+                // Two pairs that compare equal are the same two numbers, so every sort gives the same
+                // sequence and the weighted median below is summed in it on any thread count.
+                ParallelSort(lv.begin(), lv.end(), nthreads, std::less<std::pair<double, double>>());
                 double total = 0.0;
                 for (const auto &[v, w] : lv) total += w;
                 double run = 0.0;
@@ -4921,7 +4926,7 @@ RotationScaleMerge::Result RotationScaleMerge::MergeAndStats(int n_groups, bool
     std::vector<float> group_epsilon;
     std::vector<uint8_t> group_centric;
     if (wilson_test) {
-        const gemmi::GroupOps gops = x.GetSpaceGroupOrP1().operations();
+        const gemmi::GroupOps gops = Group().operations();
         group_epsilon.resize(n_groups);
         group_centric.resize(n_groups);
         ParallelChunks(n_groups, ThreadsForWork(n_groups, nthreads), [&](int lo, int hi) {
@@ -5105,7 +5110,7 @@ RotationScaleMerge::Result RotationScaleMerge::MergeAndStats(int n_groups, bool
     // --scaling-high-resolution (d_min_limit) wins; the P1 search merge (for_search) is never cut, so
     // the space-group search still sees the full range.
     const std::optional<double> effective_d_min = ApplyResolutionCutoff(
-        result.merged, d_min_limit, resolution_cutoff_method, resolution_cc_target, for_search, logger,
+        result.merged, d_min_limit, resolution_cutoff_method, resolution_cc_target, for_search, *logger,
         &result.resolution_fit_A);
 
     // The error model has to be calibrated on the reflections that are kept, not on the ones that are
@@ -5266,22 +5271,22 @@ RotationScaleMerge::Result RotationScaleMerge::MergeAndStats(int n_groups, bool
         // equivalent of and which can only ever be the more optimistic of the two.
         const auto em = ToXdsErrorModel(error_model_a, error_model_b);
         if (error_model_b_unmeasured)
-            logger.Warning("Error model (XDS convention): a={:.3f}, b NOT MEASURABLE - fewer than one "
+            logger->Warning("Error model (XDS convention): a={:.3f}, b NOT MEASURABLE - fewer than one "
                            "intensity bin's worth of reflections are strong enough to constrain it, so "
                            "b is held at 0 and ISa is not reported. chi2={:.2f}. This says the data do "
                            "not reach far enough for a systematic error to be seen, not that there is "
                            "none; a resolution range matched to the signal would measure it",
                            em.a, error_model_chi2);
         else if (!error_model_b_resolved)
-            logger.Info("Error model (XDS convention): a={:.3f} b={:.3e} - b is not resolved from zero "
+            logger->Info("Error model (XDS convention): a={:.3f} b={:.3e} - b is not resolved from zero "
                         "(under two standard errors), so ISa is undetermined (1/b = {:.1f})",
                         em.a, em.b, em.isa);
         else if (!full_stats)
             // Neither the asymptote nor the chi2 median was measured on this merge, so neither is
             // reported for it.
-            logger.Info("Error model (XDS convention): a={:.3f} b={:.3e} ISa={:.1f}", em.a, em.b, em.isa);
+            logger->Info("Error model (XDS convention): a={:.3f} b={:.3e} ISa={:.1f}", em.a, em.b, em.isa);
         else
-            logger.Info("Error model (XDS convention): a={:.3f} b={:.3e} ISa={:.1f} "
+            logger->Info("Error model (XDS convention): a={:.3f} b={:.3e} ISa={:.1f} "
                         "(strong-reflection asymptote {:.1f}) chi2={:.2f}",
                         em.a, em.b, em.isa, result.isa_asymptotic, error_model_chi2);
     }
@@ -5290,12 +5295,12 @@ RotationScaleMerge::Result RotationScaleMerge::MergeAndStats(int n_groups, bool
     // is not written needs neither. French-Wilson is deferred until after the anomalous accumulator
     // below has attached I(+)/I(-), so the two hands get their amplitudes in one pass.
     if (full_stats)
-        AssignRfreeFlags(result.merged, x.GetSpaceGroupOrP1(), rfree_fraction, 500, reference_cell);
+        AssignRfreeFlags(result.merged, Group(), rfree_fraction, 500, reference_cell);
 
     if (reject_count > 0)
-        logger.Info("Merge outlier rejection: dropped {} observations", reject_count);
+        logger->Info("Merge outlier rejection: dropped {} observations", reject_count);
     if (wilson_test && wilson.n_tested > 0)
-        logger.Info("Wilson outlier test: {} observations tested, "
+        logger->Info("Wilson outlier test: {} observations tested, "
                     "tail scale {:.2f}, bound E^2 > {:.1f} (acentric), {} rejected", wilson.n_tested,
                     wilson.tail_scale, wilson.bound, wilson.n_rejected);
     // ---- Statistics (report_shell_count shells): completeness, multiplicity, <I/sigma>, R_meas, CC1/2,
@@ -5332,7 +5337,7 @@ RotationScaleMerge::Result RotationScaleMerge::MergeAndStats(int n_groups, bool
         // than an asymmetric unit of it - so enumerating it there is the largest single piece of work in
         // the merge that nothing goes on to read.
         if (reference_cell && !for_search)
-            PossiblePerShell(x.GetSpaceGroupOrP1(), *reference_cell, grid_d_min, grid_d_max,
+            PossiblePerShell(Group(), *reference_cell, grid_d_min, grid_d_max,
                              shells, merge_friedel, possible, nthreads);
         for (int s = 0; s < n_shells; ++s) sa[s].possible = possible[s];
 
@@ -5470,8 +5475,8 @@ RotationScaleMerge::Result RotationScaleMerge::MergeAndStats(int n_groups, bool
         CorrelationCoefficient ha_cc_overall;
         bool ha_measured = false;
         if (!for_search && full_stats) {
-            const HKLKeyGenerator anom_keygen(/*merge_friedel=*/false, x.GetSpaceGroupOrP1());
-            const gemmi::GroupOps gops = x.GetSpaceGroupOrP1().operations();
+            const HKLKeyGenerator anom_keygen(/*merge_friedel=*/false, Group());
+            const gemmi::GroupOps gops = Group().operations();
             // [hand][half]: hand 0 = I(+), 1 = I(-). The half sets are what CCanom needs - the anomalous
             // difference formed twice, once from each half of the observations - and their sums add back
             // to the whole-mate I(+)/I(-) that SigAno and the export read, so there is one accumulation.
@@ -5643,7 +5648,7 @@ RotationScaleMerge::Result RotationScaleMerge::MergeAndStats(int n_groups, bool
                 // The completeness denominator counted the same way: two per acentric, one per centric.
                 if (reference_cell) {
                     std::vector<int> possible_hand(n_shells, 0);
-                    PossiblePerShell(x.GetSpaceGroupOrP1(), *reference_cell, grid_d_min, grid_d_max,
+                    PossiblePerShell(Group(), *reference_cell, grid_d_min, grid_d_max,
                                      shells, /*merge_friedel=*/false, possible_hand, nthreads);
                     for (int sh = 0; sh < n_shells; ++sh) ha[sh].possible = possible_hand[sh];
                 }
@@ -5716,7 +5721,7 @@ RotationScaleMerge::Result RotationScaleMerge::MergeAndStats(int n_groups, bool
             ho.cc_anom = overall.cc_anom;
         }
         if (own_table && std::isfinite(overall.abs_diff_over_sigma_anomalous))
-            logger.Info("Anomalous signal SigAno = {:.2f} ({} acentric pairs), CCanom = {:.3f} ({} pairs "
+            logger->Info("Anomalous signal SigAno = {:.2f} ({} acentric pairs), CCanom = {:.3f} ({} pairs "
                         "split in both hands)",
                         overall.abs_diff_over_sigma_anomalous, sig_n, overall.cc_anom, cc_anom_n);
         return out;
@@ -5798,7 +5803,7 @@ RotationScaleMerge::Result RotationScaleMerge::MergeAndStats(int n_groups, bool
     // Attach the per-reflection anomalous split so the writer can emit I(+)/I(-) by default (each merged
     // reflection maps to its Friedel-ASU key; in an anomalous merge both mates map to the same key).
     if (!anom_export.empty()) {
-        const HKLKeyGenerator anom_keygen(/*merge_friedel=*/false, x.GetSpaceGroupOrP1());
+        const HKLKeyGenerator anom_keygen(/*merge_friedel=*/false, Group());
         for (auto &r : result.merged) {
             const HKLKey ak = anom_keygen(r.h, r.k, r.l);
             const auto it = anom_export.find(HKLKey{ak.h, ak.k, ak.l, true}.pack());
@@ -5812,10 +5817,10 @@ RotationScaleMerge::Result RotationScaleMerge::MergeAndStats(int n_groups, bool
     if (full_stats) {
         FrenchWilsonOptions fw_opts;
         fw_opts.num_threads = static_cast<int>(nthreads);
-        ApplyFrenchWilson(result.merged, x.GetSpaceGroupOrP1(), fw_opts);
+        ApplyFrenchWilson(result.merged, Group(), fw_opts);
     }
 
-    logger.Info("Merge complete ({} unique reflections)", result.merged.size());
+    logger->Info("Merge complete ({} unique reflections)", result.merged.size());
     return result;
 }
 
@@ -5829,10 +5834,14 @@ namespace {
     // rocking events per frame the partials' scale read ISa 10.5 against 8.7 from the fulls alone, and
     // SHELXL R1 0.105 against 0.062 - the partiality-model error the partial scale takes up is shared by
     // symmetry mates measured at the same rocking geometry, so their agreement cannot see it.
+    bool ErrorModelResolved(const RotationScaleMerge::Result &r) {
+        return r.isa_resolved && std::isfinite(r.isa) && r.isa > 0.0;
+    }
+
     bool PartialScalingWins(const RotationScaleMerge::Result &with, const RotationScaleMerge::Result &without,
                             bool prior) {
-        const bool w_ok = with.isa_resolved && std::isfinite(with.isa) && with.isa > 0.0;
-        const bool o_ok = without.isa_resolved && std::isfinite(without.isa) && without.isa > 0.0;
+        const bool w_ok = ErrorModelResolved(with);
+        const bool o_ok = ErrorModelResolved(without);
         const bool prior_ok = prior ? w_ok : o_ok;
         const bool other_ok = prior ? o_ok : w_ok;
         return prior_ok || !other_ok ? prior : !prior;
@@ -5840,7 +5849,10 @@ namespace {
 }   // namespace
 
 RotationScaleMerge::Result RotationScaleMerge::Run(bool for_search, bool full_stats,
-                                                  bool measure_cc_before_corrections) {
+                                                  bool measure_cc_before_corrections,
+                                                  const gemmi::SpaceGroup *group) {
+    run_group = group;
+    struct ClearGroup { const gemmi::SpaceGroup *&g; ~ClearGroup() { g = nullptr; } } clear_group{run_group};
     // Whether a frame's scale can be fitted on its partials depends on whether the partials tell a change
     // of scale from an error of the partiality model. The count of rocking events per frame is the prior
     // (MIN_EVENTS_PER_FRAME_SCALE_PARTIALS), but the counts overlap between small molecules (2.6-43) and
@@ -5850,17 +5862,33 @@ RotationScaleMerge::Result RotationScaleMerge::Run(bool for_search, bool full_st
     // where the other did not (PartialScalingWins); later merges keep that choice. A search merge is in
     // P1 or a subgroup, where a reflection has two or three observations, so it takes the prior. The arm
     // returned is the one run last, so that every per-pass state the merge leaves behind is that arm's.
+    //
+    // Where the prior is to scale the partials and that merge resolves its error model, the other arm
+    // cannot change the choice, so it is not made. The partials' arm then starts from G = 1, which is
+    // what the fulls-alone arm would have left it (it is the only state that arm leaves for the next),
+    // so the merge is the one the pair would have returned.
     if (rocking_event_frames_at_start < 0)
         rocking_event_frames_at_start = RockingEventFrames(&rocking_events_at_start);
     const double events_per_frame = n_frames > 0 ? static_cast<double>(rocking_events_at_start) / n_frames : 0.0;
     const bool prior = events_per_frame >= MIN_EVENTS_PER_FRAME_SCALE_PARTIALS;
     if (for_search && !scale_partials_choice)
         return RunArm(prior, for_search, full_stats, measure_cc_before_corrections);
+    if (!scale_partials_choice && prior) {
+        std::fill(g_partial.begin(), g_partial.end(), 1.0);
+        Result with = RunArm(true, for_search, full_stats, measure_cc_before_corrections);
+        if (ErrorModelResolved(with)) {
+            scale_partials_choice = true;
+            logger->Info("Per-frame scale: fitted on the partials ({:.1f} rocking events per frame) - ISa {:.1f} "
+                        "with the partials scaled; the fulls-alone merge is made only when this one has no "
+                        "resolved error model", events_per_frame, with.isa);
+            return with;
+        }
+    }
     if (!scale_partials_choice) {
         Result without = RunArm(false, for_search, full_stats, measure_cc_before_corrections);
         Result with = RunArm(true, for_search, full_stats, measure_cc_before_corrections);
         scale_partials_choice = PartialScalingWins(with, without, prior);
-        logger.Info("Per-frame scale: {} ({:.1f} rocking events per frame) - ISa {:.1f}{} with the partials "
+        logger->Info("Per-frame scale: {} ({:.1f} rocking events per frame) - ISa {:.1f}{} with the partials "
                     "scaled, {:.1f}{} from the fulls alone",
                     *scale_partials_choice ? "fitted on the partials" : "from the fulls alone",
                     events_per_frame, with.isa, with.isa_resolved ? "" : " (unresolved)",
@@ -5873,7 +5901,7 @@ RotationScaleMerge::Result RotationScaleMerge::Run(bool for_search, bool full_st
 
 RotationScaleMerge::Result RotationScaleMerge::RunArm(bool scale_partials, bool for_search, bool full_stats,
                                                      bool measure_cc_before_corrections) {
-    HKLKeyGenerator keygen(merge_friedel, x.GetSpaceGroupOrP1());
+    HKLKeyGenerator keygen(merge_friedel, Group());
 
     // Start from the corr Ingest built, so this pass runs the scaling_iter iterations it was asked for
     // rather than continuing the previous pass's (see the header).
@@ -5918,7 +5946,7 @@ RotationScaleMerge::Result RotationScaleMerge::RunArm(bool scale_partials, bool
         smooth_window = std::max(static_cast<int>(std::lround(smooth_g_deg / osc_deg)),
                                  SMOOTH_G_MIN_ROCKING_EVENTS * rocking_event_frames_at_start);
         if (smooth_window % 2 == 0) ++smooth_window;
-        logger.Info("Per-frame scale smoothing window: {} frames ({:.1f} deg; {} deg asked, rocking "
+        logger->Info("Per-frame scale smoothing window: {} frames ({:.1f} deg; {} deg asked, rocking "
                     "curves span {} frames)", smooth_window, smooth_window * osc_deg, smooth_g_deg,
                     rocking_event_frames_at_start);
     }
@@ -5952,7 +5980,7 @@ RotationScaleMerge::Result RotationScaleMerge::RunArm(bool scale_partials, bool
         frame_scaled_scratch.assign(n_frames, 1);
         partial_loop.converged = true;
     } else if (!scaled_on_gpu) {
-        const PartialLoopKey key{x.GetSpaceGroupOrP1().xhm(), merge_friedel, d_min_limit, d_max_limit,
+        const PartialLoopKey key{Group().xhm(), merge_friedel, d_min_limit, d_max_limit,
                                  min_partiality, smooth_window, scaling_iter};
         const auto memo = std::find_if(partial_loop_memos.begin(), partial_loop_memos.end(),
                                        [&](const PartialLoopMemo &m) { return m.key == key; });
@@ -5966,7 +5994,7 @@ RotationScaleMerge::Result RotationScaleMerge::RunArm(bool scale_partials, bool
                            [&](int lo, int hi) {
                 for (int i = lo; i < hi; ++i) partials[i].corr = m.corr[i];
             });
-            logger.Info("Per-frame scaling (partials): the same grouping and settings as an earlier merge, "
+            logger->Info("Per-frame scaling (partials): the same grouping and settings as an earlier merge, "
                         "so its {} iterations stand", partial_loop.iterations);
         } else {
             // The loop works on its own compact copy of the partials (ScalingObs in the header) and hands
@@ -6125,7 +6153,7 @@ RotationScaleMerge::Result RotationScaleMerge::RunArm(bool scale_partials, bool
             n_dropped = dropped;
         }
         if (n_dropped > 0)
-            logger.Info("Space-group search: ignoring {} observations with |zeta| < {:.2f} "
+            logger->Info("Space-group search: ignoring {} observations with |zeta| < {:.2f} "
                         "(they cross the Ewald sphere near-tangentially and are measured worst)",
                         n_dropped, search_min_zeta);
     }
@@ -6162,7 +6190,7 @@ RotationScaleMerge::Result RotationScaleMerge::RunArm(bool scale_partials, bool
                         for (int i = lo; i < hi; ++i)
                             if (reject[partials[i].frame]) partials[i].corr = 0.0f;
                     });
-                logger.Info("Space-group search: ignoring {} of {} frames whose scale came out below "
+                logger->Info("Space-group search: ignoring {} of {} frames whose scale came out below "
                             "1/{:.0f} of the run's typical frame (the crystal barely diffracted there)",
                             n_rejected, n_frames, 1.0 / SEARCH_MIN_SCALE_RATIO);
             }
@@ -6194,7 +6222,7 @@ RotationScaleMerge::Result RotationScaleMerge::RunArm(bool scale_partials, bool
                     for (int i = lo; i < hi; ++i)
                         if (reject[partials[i].frame]) partials[i].corr = 0.0f;
                 });
-            logger.Info("Rejected {} of {} frames correlating below {:.2f} with the merged reference",
+            logger->Info("Rejected {} of {} frames correlating below {:.2f} with the merged reference",
                         n_rejected, n_frames, min_cc_for_image);
         }
     }
@@ -6279,7 +6307,7 @@ RotationScaleMerge::Result RotationScaleMerge::RunArm(bool scale_partials, bool
                 o.zeta = 0.0f; o.delta_phi = 0.0f; o.bkg = 0.0f;
             }
         });
-        logger.Info("3D combine{} (GPU): {} fulls, {} rocking events dropped for a saturated pixel",
+        logger->Info("3D combine{} (GPU): {} fulls, {} rocking events dropped for a saturated pixel",
                     scaled_fulls_on_gpu ? " + scale-fulls" : "", nf, fulls_dropped_overloaded);
         combined_on_gpu = true;
     }
@@ -6309,7 +6337,7 @@ RotationScaleMerge::Result RotationScaleMerge::RunArm(bool scale_partials, bool
                 if (apply[o.frame] && std::isfinite(o.corr))
                     o.corr = static_cast<float>(o.corr * ratio[o.frame]);
         });
-        logger.Info("Scaled fulls (XDS order, Unity model)");
+        logger->Info("Scaled fulls (XDS order, Unity model)");
     }
     // The total scale each frame's fulls were amplified by - partial scale, incident flux and the fulls'
     // own G (corr = 1/G here) - taken BEFORE the correction surfaces fold per-observation terms into
@@ -6373,7 +6401,7 @@ RotationScaleMerge::Result RotationScaleMerge::RunArm(bool scale_partials, bool
         cc_half_before_corrections = MergeAndStats(n_groups, for_search,
                                                    combined_on_gpu && scaled_fulls_on_gpu,
                                                    /*full_stats=*/false).statistics.overall.cc_half;
-        logger.Info("CC1/2 before the correction surfaces: {:.4f} (the merge just above; the pass this "
+        logger->Info("CC1/2 before the correction surfaces: {:.4f} (the merge just above; the pass this "
                     "one is judged against fits no surfaces, so this is what the two are compared on)",
                     *cc_half_before_corrections);
     }

image_analysis/scale_merge/RotationScaleMerge.h

@@ -132,7 +132,11 @@ public:
     //   data needs it (the pre-pass fits no surfaces, so only the uncorrected number is the same
     //   measurement on both sides); it is a whole extra merge, so the offline re-scale path, which
     //   compares nothing, asks for it to be left out.
-    Result Run(bool for_search, bool full_stats, bool measure_cc_before_corrections);
+    // group: the space group to merge in, for this Run only; null takes the experiment's. Given, the merge
+    //   reads the experiment for nothing that changes during a run, so a caller may merge in its own group
+    //   (the P1 cross-check) on another thread while it goes on using the experiment.
+    Result Run(bool for_search, bool full_stats, bool measure_cc_before_corrections,
+               const gemmi::SpaceGroup *group = nullptr);
 
     // Whether Run() writes the per-frame G / CC / mosaicity back onto the outcomes (on by default). Off
     // for a merge that is not the run's answer - the P1 cross-check - so the per-image table and the
@@ -143,6 +147,10 @@ public:
     // it is a copy as long as the fulls, wanted only by a caller that writes them.
     void SetExportScaledFulls(bool on) { export_scaled_fulls = on; }
 
+    // Where the engine logs from now on (by default the logger it was made with) - a merge run beside
+    // other work logs into a held buffer, so the run's log keeps one order.
+    void SetLogger(Logger &to) { logger = &to; }
+
     // Override the high-resolution cut for the next Run() - used to gate the de-novo P1 search pass at
     // <I/sigma> >= 1 without cutting the final in-symmetry merge. Reset to the manual limit afterwards.
     void SetDMinLimit(std::optional<double> d_min_A) { d_min_limit = d_min_A; }
@@ -201,8 +209,11 @@ private:
     std::vector<IntegrationOutcome> &partials_out;   // written back at the end of scaling
     std::optional<UnitCell> reference_cell;
     size_t nthreads;
-    Logger &logger;
+    Logger *logger;
     std::string observation_dump_path;
+    // The group the current Run merges in (see Run); null outside a Run that was given one.
+    const gemmi::SpaceGroup *run_group = nullptr;
+    const gemmi::SpaceGroup &Group() const { return run_group ? *run_group : x.GetSpaceGroupOrP1(); }
 
     // Fixed settings snapshot (read once in the ctor).
     int n_frames = 0;

rugnux/Rugnux.cpp

@@ -8854,6 +8854,38 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b
         if (superseded)
             logger.Info("This pass is going to be re-run, so it makes no report and writes no files.");
 
+        // Whether the P1 cross-check dataset is made (see where it is written, below).
+        const bool p1_integration_complete =
+                search_space_group || (indexer && indexer->GetPredictionCentring() == 'P');
+        const bool p1_crosscheck = result.consensus_cell && write_files && config_.write_merged
+                                   && !geometry_prepass && !superseded && config_.write_p1_crosscheck
+                                   && p1_integration_complete && is_rotation;
+        // The P1 cross-check merge is made from here on, beside everything that comes before it is
+        // written. Nothing here reads what it writes or writes what it reads: the engine merges its own
+        // copy of the observations, in a group of its own rather than the experiment's, and writes no
+        // per-frame scale back onto the outcomes; and nothing up to take_p1 calls the engine, changes
+        // the outcomes it was ingested from, or changes the experiment. Its log lines are held and
+        // printed where it is taken.
+        std::optional<RotationScaleMerge::Result> p1_merged_early;
+        Logger p1_log = Logger::Buffered();
+        std::future<RotationScaleMerge::Result> p1_ahead;
+        if (p1_crosscheck && rsm) {
+            rsm->SetWriteBackPerFrameScale(false);
+            rsm->SetExportScaledFulls(false);
+            rsm->SetLogger(p1_log);
+            p1_ahead = std::async(std::launch::async, [&rsm] {
+                return rsm->Run(/*for_search=*/false, /*full_stats=*/true,
+                                /*measure_cc_before_corrections=*/false, &gemmi::get_spacegroup_p1());
+            });
+        }
+        auto take_p1 = [&] {
+            if (!p1_ahead.valid()) return;
+            p1_merged_early = p1_ahead.get();
+            rsm->SetLogger(logger);
+            rsm->SetWriteBackPerFrameScale(true);
+            p1_log.ReplayInto(logger);
+        };
+
         const auto &twin_sg_opt = experiment_.GetGemmiSpaceGroup();
         const gemmi::SpaceGroup *twin_sg = twin_sg_opt ? &*twin_sg_opt : nullptr;
 
@@ -9307,14 +9339,6 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b
             }
         }
 
-        // Whether the P1 cross-check dataset is made (see where it is written, below).
-        const bool p1_integration_complete =
-                search_space_group || (indexer && indexer->GetPredictionCentring() == 'P');
-        const bool p1_crosscheck = result.consensus_cell && write_files && config_.write_merged
-                                   && !geometry_prepass && !superseded && config_.write_p1_crosscheck
-                                   && p1_integration_complete && is_rotation;
-        std::optional<RotationScaleMerge::Result> p1_merged_early;
-
         // Model validation runs BEFORE the reflection files are written, because it is what settles the
         // frame they are written in: the enantiomorph, which merged intensities cannot choose, and -
         // with no reference MTZ - the alternative indexing. Both are relabelings of the same
@@ -9339,28 +9363,11 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b
                                             static_cast<size_t>(config_.nthreads),
                                             wavelength, report_shell_d_min);
             };
-            // The P1 cross-check merge and this validation read nothing the other writes, so they
-            // are made at the same time. The validation reads the merge above, the cell and the group,
-            // none of which the scaling engine touches; the engine merges its own copy of the
-            // observations, in P1, and writes back only per-image scale, CC and mosaicity, which
-            // nothing before the cross-check reads (the unmerged MTZ replaces the mosaicity it
-            // carries after the merge anyway). Every relabelling the validation decides reaches the
-            // P1 merge afterwards, through merge_to_written, as it did when it was made below. The
-            // validation's log lines are held and printed as one block once it is done.
-            ModelValidationResult validation;
-            if (p1_crosscheck && rsm) {
-                Logger held = Logger::Buffered();
-                auto pending = std::async(std::launch::async, validate, std::ref(held));
-                experiment_.SpaceGroupNumber(1);
-                rsm->SetWriteBackPerFrameScale(false);
-                p1_merged_early = rsm->Run(/*for_search=*/false, /*full_stats=*/true,
-                                           /*measure_cc_before_corrections=*/false);
-                rsm->SetWriteBackPerFrameScale(true);
-                experiment_.SetSpaceGroup(data_sg);
-                validation = pending.get();
-                held.ReplayInto(logger);
-            } else
-                validation = validate(logger);
+            // The P1 cross-check merge (started above) is taken before anything below acts on what the
+            // validation decides: relabelling the outcomes and the group is what it must not see happen.
+            // Every such relabelling reaches the P1 merge afterwards, through merge_to_written.
+            ModelValidationResult validation = validate(logger);
+            take_p1();
             // A model that was asked for and could not be used has to say so where anyone will see
             // it. Without this the run ends successfully with no R-free, no maps and nothing in the
             // report - indistinguishable from a run that was never given --model at all.
@@ -9556,6 +9563,7 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b
                 }
                 // Both the merge and the MTZ read the group from the experiment, so it is set for
                 // the whole of it and restored after.
+                take_p1();
                 experiment_.SpaceGroupNumber(1);
                 if (!p1_merged_early)
                     rsm->SetWriteBackPerFrameScale(false);