mx/Jungfraujoch
3
2
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases 81 Wiki Activity

Two-pass: a partial pre-pass is a sweep of the frames it read, and the guard compares those frames

Browse Source 
With --prepass-fraction below 1 the geometry pre-pass stopped integrating at the fraction but still
handed the merge and the post-refinement an outcome vector the length of the whole sweep. The
fulls' per-frame log-scale smoother (Whittaker, lambda by cross-validation) then carried its
curve straight across the empty tail; on a sparse small-molecule sweep with 900 empty frames the
extrapolated scale fell further each round until it underflowed to zero, log G became -inf there,
the cross-validation error became 0*inf = NaN for every lambda, lambda fell back to 0, the
smoother's band matrix was singular on the frames without fulls, and every per-frame scale in the
sweep went NaN. The merge was then empty and the run aborted with "RotationScaleMerge: resolution
calculation failed". The outcome vector is now cut to the frames the pre-pass read, so its merge,
its rocking-events-per-frame prior, its smoothed mosaicity and its post-refinement all describe a
sweep of those frames. The pre-pass mosaicity that the second pass predicts with is carried to the
frames after the pre-pass at the median of what the pre-pass fitted (on a small-molecule sweep,
predicting the second half at each frame's own estimate alone lost a fifth of the predictions and
took R_meas from 5.2 % to 7.2 %).

The two-pass quality guard compares the refined pass with the pre-pass on their search merges.
Against a pre-pass of the first half, a merge of the whole sweep always holds more reflections,
more signal and more of the axial rows, so the guard could never fire (it fired on 6 battery sets
at fraction 1 and on none at 0.5). A refined pass judged against a partial pre-pass now fills the
guard's numbers from its own integration of the same frames, merged the way the pre-pass merged
them (P1 or the fixed group, search merge, no surfaces). At fraction 1 nothing changes.

Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01SVmAWnzCmRKAXVUCdc4iNi
This commit is contained in:
leonarski_f
2026-10-07 12:16:09 +02:00
co-authored by Claude Opus 5.5
parent ed2e0ffadf
commit 01a4e95041
1 changed files with 71 additions and 13 deletions
Download Patch File Download Diff File Expand all files Collapse all files
rugnux/Rugnux.cpp
+71 -13
View File
@@ -6281,6 +6281,12 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b
}
};
// The geometry pre-pass integrates only the first prepass_fraction of the sweep (pass 2 reads it all).
const int prepass_end = start_image
+ static_cast<int>(std::lround((end_image - start_image) * config_.prepass_fraction));
const int prepass_images = (prepass_end - start_image + config_.stride - 1) / config_.stride;
const int loop_end = geometry_prepass ? prepass_end : end_image;
auto full_worker = [&]() {
pin_gpu(); // round-robin per worker thread; must precede engine construction
// Before the analysis engine, which page-locks these bytes for its uploads and unregisters
@@ -6291,9 +6297,6 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b
AzimuthalIntegrationProfile profile(mapping);
HarmonicEvidence worker_harmonic;
std::vector<FootprintOffset> worker_offsets;
const int loop_end = geometry_prepass
? start_image + static_cast<int>(std::lround((end_image - start_image) * config_.prepass_fraction))
: end_image;
while (!cancelled_) {
const int ordinal = next_ordinal.fetch_add(1);
@@ -6424,6 +6427,14 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b
Note(fmt::format("{} images, {:.0f} per second", result.images_processed,
result.images_processed / std::max(1e-6, result.image_loop_time_s)));
// A pre-pass that stopped at loop_end is a sweep of the frames it read, and that is all the merge and
// the post-refinement are handed. Left in, the empty tail is a stretch of frames every per-frame
// smoothing in the merge extrapolates over: measured, the fulls' log-scale curve, carried straight
// across 900 empty frames, fell further every round until the scale there was zero, and the
// cross-validation of its smoothness then took every scale in the sweep to NaN.
if (geometry_prepass)
indexer->GetIntegrationOutcome().resize(prepass_images);
// The footprint of the reflection rather than of the spot: the pre-scan's widths plus where the
// spots sit against their predictions (SpotFootprint.h). Like the widths it acts only where the
// reflection outgrows the r1 disk, and like them it is held back for the canonical pass.
@@ -6984,14 +6995,27 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b
const auto &rot_ss = experiment_.GetScalingSettings();
const bool is_rotation = experiment_.IsRotationIndexing(); // rotation indexing -> rotation scaling/merge
std::optional<RotationScaleMerge> rsm;
std::optional<ScaleMergeResult> guard_merge; // what the quality guard reads, where not the search merge
std::optional<PostRefineObservations> prepass_postrefine_obs;
// The rotation geometry post-refinement (see its call sites below for what it is for).
const auto post_refine_geometry = [&] {
const auto &outcomes = indexer->GetIntegrationOutcome();
if (geometry_prepass) {
prepass_mosaicity_.assign(outcomes.size(), NAN);
std::vector<float> fitted;
for (size_t o = 0; o < outcomes.size(); ++o)
if (outcomes[o].mosaicity_deg) prepass_mosaicity_[o] = *outcomes[o].mosaicity_deg;
if (outcomes[o].mosaicity_deg) {
prepass_mosaicity_[o] = *outcomes[o].mosaicity_deg;
fitted.push_back(*outcomes[o].mosaicity_deg);
}
// A pre-pass of part of the sweep has fitted nothing for the frames after it, and they are
// predicted at the median of what it did fit, as the merge fills its own frames without a
// fit. Left at the frame's own estimate alone, the second half of a small-molecule sweep
// predicted a fifth fewer reflections and merged at R_meas 7.2 % against 5.2 %.
if (!fitted.empty() && outcomes.size() < static_cast<size_t>(images_to_process)) {
std::nth_element(fitted.begin(), fitted.begin() + fitted.size() / 2, fitted.end());
prepass_mosaicity_.resize(images_to_process, fitted[fitted.size() / 2]);
}
}
if (result.consensus_cell.has_value()) {
if (const auto rot = indexer->FinalizeRotationIndexing(); rot && rot->axis) {
@@ -7118,6 +7142,38 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
"Rotation scaling/merging (RotationScaleMerge) does not support "
"wedge refinement");
// The quality guard judges a refined pass against the pre-pass on the search merge, and a
// pre-pass that integrated only the first part of the sweep merged only that part. Against
// it a merge of the whole sweep holds more reflections, more signal and more of the axial
// rows for no reason but the frames, and the guard could never find it worse. So the guard
// reads this pass's own integration of the same frames, merged the way the pre-pass merged
// them. The merge writes its per-frame scale back, so it gets copies of the outcomes and
// only borrows the reflections.
if (quality_guard_pass1_ && prepass_images < images_to_process) {
logger.Info("Two-pass: the quality guard reads this pass's merge of the first {} images, "
"the ones the pre-pass integrated", prepass_images);
auto &outcomes = indexer->GetIntegrationOutcome();
std::vector<IntegrationOutcome> first_part(prepass_images);
for (int i = 0; i < prepass_images; ++i) {
ReflectionVector reflections = std::move(outcomes[i].reflections);
first_part[i] = outcomes[i];
first_part[i].reflections = std::move(reflections);
}
{
RotationScaleMerge first_part_rsm(experiment_, first_part, result.consensus_cell,
static_cast<int>(config_.scaling_iter), config_.nthreads,
logger, "");
first_part_rsm.Ingest();
first_part_rsm.SetExportScaledFulls(false);
first_part_rsm.SetFrenchWilson(false);
auto r = first_part_rsm.Run(!experiment_.GetGemmiSpaceGroup().has_value(),
/*full_stats=*/false, /*measure_cc_before_corrections=*/true);
guard_merge = ScaleMergeResult{std::move(r.merged), std::move(r.statistics),
r.cc_half_before_corrections};
}
for (int i = 0; i < prepass_images; ++i)
outcomes[i].reflections = std::move(first_part[i].reflections);
}
// A reference MTZ is allowed for rotation: it fixes the space group / cell (on the CLI) and
// resolves the indexing ambiguity (below), but is NOT used to scale - the rotation merge stays
// self-consistent.
@@ -7279,15 +7335,17 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b
// The two-pass quality guard judges the second pass against the first on THIS merge, which both
// passes make in the same terms - P1 (or the group the user fixed), full range, no correction
// surfaces - where their final merges need not even be in the same group. See RunAllPasses.
// surfaces - where their final merges need not even be in the same group. See RunAllPasses. Where
// the pre-pass integrated only part of the sweep, it is the same merge of those frames (above).
{
const auto &o = sm.statistics.overall;
const ScaleMergeResult &counted = guard_merge ? *guard_merge : sm;
const auto &o = counted.statistics.overall;
result.search_merge_completeness_measured = o.possible_unique_reflections > 0;
result.search_merge_completeness = result.search_merge_completeness_measured
? 100.0 * o.unique_reflections / o.possible_unique_reflections : 0.0;
result.search_merge_cc_half = sm.cc_half_before_corrections;
if (!sm.statistics.shells.empty())
result.search_merge_inner_cc_half = sm.statistics.shells.front().cc_half;
result.search_merge_cc_half = counted.cc_half_before_corrections;
if (!counted.statistics.shells.empty())
result.search_merge_inner_cc_half = counted.statistics.shells.front().cc_half;
// ... and the reflections it holds on the principal axial rows, which the guard needs
// separately because CC1/2 cannot express them (see RunAllPasses). Only the low-order
// part is counted: a screw is read off the strong start of its row - the evidence gate
@@ -7316,7 +7374,7 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b
};
constexpr float AXIAL_ROW_LOW_ORDER_D_A = 4.0f;
result.search_merge_axial_reflections = std::count_if(
sm.merged.begin(), sm.merged.end(), [&](const MergedReflection &r) {
counted.merged.begin(), counted.merged.end(), [&](const MergedReflection &r) {
return ((r.h == 0) + (r.k == 0) + (r.l == 0)) == 2
&& r.d > AXIAL_ROW_LOW_ORDER_D_A && centring_allows(r);
});
@@ -7324,14 +7382,14 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b
// what the guard decides on.
constexpr float STRONG_I_OVER_SIGMA = 2.0f;
result.search_merge_strong_reflections = std::count_if(
sm.merged.begin(), sm.merged.end(), [&](const MergedReflection &r) {
counted.merged.begin(), counted.merged.end(), [&](const MergedReflection &r) {
return r.sigma > 0 && r.I / r.sigma >= STRONG_I_OVER_SIGMA && centring_allows(r);
});
result.search_merge_negative_reflections = std::count_if(
sm.merged.begin(), sm.merged.end(), [&](const MergedReflection &r) {
counted.merged.begin(), counted.merged.end(), [&](const MergedReflection &r) {
return r.sigma > 0 && r.I / r.sigma <= -STRONG_I_OVER_SIGMA && centring_allows(r);
});
result.search_merge_reflections = std::count_if(sm.merged.begin(), sm.merged.end(),
result.search_merge_reflections = std::count_if(counted.merged.begin(), counted.merged.end(),
centring_allows);
}