From 8153646716b9d452db0a56eeb8599d1ed10b1513 Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Fri, 25 Sep 2026 08:41:31 +0200 Subject: [PATCH 1/4] rugnux: geometry walk keeps a round the validation spots prefer, and probes small moves The between-pass geometry walk kept a round only when the realised held-out residual fell by more than its noise, and started only on a fit move of a trust-region step or more. The residual is dominated by low-resolution reflections, where a distance and the compensating cell scale move every spot alike, and its centroids are taken inside a disc centred on the prediction, so it barely sees a distance error that costs the high-resolution shells. On 8pqd the canonical pass ran at 96.456 mm; the fit asked for 95.878 mm (0.6 %, less than a step, so no walk). Forced, that round read the residual only 0.74 sigma lower but put 70.4 % of the validation spots on the lattice against 58.1 %. Now a round is kept when either the residual falls beyond its noise (HeldOutResidualFell) or the validation evidence prefers it (ValidationEvidencePrefers, z = 3.29). A move of less than a step is first tried as two index-only probes on the validation frames (fit's geometry vs the one in hand) and pays for a re-integrated round only where the fit's geometry scores higher; walks started by a large move run as before. 8pqd: 96.456 -> 95.878 mm, d_min 1.374 -> 1.306 A, R_free .218 -> .206, REFMAC R_free .209 -> .202, Wilson B 34.4 -> 30.7 (pool had 1.325 A / .209 / .203). 6vww 7n0i 7ris 8egn 8xtf 9qw8 9w3y 6qaj 6cdl 9i0a myob_x06da_split lyso_x06da_half_image unchanged (probes say the geometry in hand stands); probe cost 1-10 s per run. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01D1G8gJVAy6gp1K5Dz3NE5C --- docs/CHANGELOG.md | 1 + docs/CPU_DATA_ANALYSIS_INDEXING.md | 2 + rugnux/Rugnux.cpp | 100 +++++++++++++++++++++++++---- 3 files changed, 90 insertions(+), 13 deletions(-) diff --git a/docs/CHANGELOG.md b/docs/CHANGELOG.md index 1ad3e1b3e..e4bbd0687 100644 --- a/docs/CHANGELOG.md +++ b/docs/CHANGELOG.md @@ -9,6 +9,7 @@ * Rugnux accepts a detector-modulation or absorption correction surface when held-out data support it on Fisher's z, and fits them in the order modulation, time, goniometer frame; this recovers corrections that were wrongly refused. * Rugnux settles a disagreement between the file's and the measured beam centre over an axis harmonic in either direction, so a file centre off along the spindle no longer leaves the run on a doubled axis. * Rugnux commits a post-refined geometry only when the held-out residual falls by more than its own noise, instead of by a fixed 2 %. +* Rugnux keeps refining the detector geometry between passes where the fit's new geometry puts more validation spots on the lattice, not only where the held-out residual falls. * Rugnux also reports `R_MEAS_WEIGHTED`, the R_meas with each observation weighted as the merge weights it, so weak frames kept in the merge at low weight can be told from real disagreement; `R_MEAS` is unchanged. * Rugnux no longer stops per-frame scaling as unsettled because of a single frame that is about to be dropped as blank. * `rugnux --polarization` is documented as the polarization degree (XDS `FRACTION_OF_POLARIZATION` = (1 + p)/2); the default 0.99 suits undulators. diff --git a/docs/CPU_DATA_ANALYSIS_INDEXING.md b/docs/CPU_DATA_ANALYSIS_INDEXING.md index 4ac0b252d..aa5b55c76 100644 --- a/docs/CPU_DATA_ANALYSIS_INDEXING.md +++ b/docs/CPU_DATA_ANALYSIS_INDEXING.md @@ -331,6 +331,8 @@ The refinement above (§7.2) runs per image against that image's spots. For rota **Whether the data determine the distance at all is asked, not assumed.** At a detector far enough away that no reflection reaches more than a few degrees of $2\theta$, a longer distance and a larger cell move every spot the same way to first order — the difference is of order $\sin^2\theta$ of the spot's own position, about 0.4 px rms per per cent of distance over a 2M detector at 820 mm against 2–3 px at the distances a crystal is usually collected at. The joint fit then finds a distance/cell pair that fits its own spot positions a little better than the header, commits it, and the pass re-integrated there finds the next pair: a walk along the degenerate direction that the realised residual never ratifies (measured on such a sweep: a header at 820 mm walked to 846 mm with the cell 3.3 % too large, the realised held-out residual flat at every round). So the same fit is asked once more with the distance **held at the header**, every other block as free as before — the nested hypothesis "the header distance is right" — and the two are compared on the one residual family that can tell them apart: the **excitation** residual. It never involves the detector, so it is blind to the distance itself; what it sees is the cell scale, and a held fit at a wrong header distance is forced into a wrong cell scale by the spot positions, which the rocking angles then refuse (measured: a header 1.4 % long leaves the held fit's excitation residual seventeen times the free fit's). Where freeing the distance lowers the held-out excitation residual below the held fit's by more than that residual's own standard error, the free fit is committed exactly as before; where it does not, the held fit is — header distance, refined beam, cell, orientation and axis — and the report says so (`POSTREFINE_DISTANCE_HELD`). The positional residual is deliberately not consulted for this: it is the family whose in-fit gain along the degenerate direction re-integration erases, and pooled with the excitation family it either drowns a decisive excitation gain in its own noise (a 54 % excitation gain read as 9 % pooled against a 9 % noise) or lends the degenerate direction a gain that is not there. Nothing is tuned here: the only input is the standard error of the residual itself, the same noise the geometry walk's rounds have to beat. The wavelength is never refined on a single crystal for the same reason in its exact form: it scales the spot positions and the rocking angles identically to the cell, so no sweep can tell the two apart at any $2\theta$. 3. **Pass 2** re-indexes de novo and re-integrates at the committed geometry. Only the **detector distance and beam centre** carry over: the refined cell, orientation and axis are what make the distance identifiable, but pass 2 re-indexes from scratch, so they are not propagated. Where that re-index indexes too few frames the run falls back to pass 1's lattice and integrates it at the refined geometry — and the cell is then **scaled to the distance it will be used at**, since a real-space cell is measured against the distance its spots were seen at, and carrying it across a distance change otherwise scales the whole cell by the ratio of the two. The orientation is untouched. + Pass 2 measures the post-refinement again, and where it still moves the geometry the run **walks**: it re-indexes and re-integrates at what the fit asks for, and repeats. A round is kept only for what it *realises*, not for what the fit predicts, and it can realise a gain in two ways, either of which has to beat its own noise: a lower held-out residual (the standard errors of the two means combined), or a larger share of the validation spots on the lattice beyond chance (the binomial noise of the two shares, z = 3.29). The residual alone misses exactly the errors that cost resolution: it is dominated by the low-resolution reflections, where a distance and the cell scale that compensates it move every spot alike, and its centroids are taken inside a disc centred on the prediction, so they follow the prediction part of the way; the high-resolution validation spots are the first to leave the lattice (measured: 0.6 % of distance read 0.74 of the residual's noise and 70.4 % against 58.1 % of the validation spots, and cost 0.07 Å of resolution). A move of a trust-region step or more starts the walk outright. A smaller one is first tried as two indexing probes on the validation frames — at the fit's geometry and at the one in hand, each stopping once the lattice is scored — and pays for a re-integrated round only where the fit's geometry scores higher. The run keeps the best round it reached. + The space group is determined **after** pass 2, on the geometry the run refined, and pass 1 does not search at all: a decision taken on the worse of the two passes and then carried forward is a constraint on the better one, and would have to be reconciled with what pass 2 later found. The guard that chooses which pass is written compares each pass's **first** merge — $P1$ on both sides, full resolution range, before the correction surfaces — which both passes produce anyway, so it never compares statistics computed in two different space groups. What it compares there is the **signal each pass measured**: the count of unique reflections merged at $I/\sigma \ge 2$. Self-consistency cannot do this job — against an external arbiter the signal count named the more accurate geometry on 23 of 27 arm-dataset pairs where $R_\mathrm{meas}$, $CC_{1/2}$ and ISa managed 13, a coin flip — and that merge's own $CC_{1/2}$ least of all, being pooled over the whole range, uncut and uncorrected, so the shells with no signal in them dominate it and they are exactly the shells a geometry move disturbs (it reads 0.13 on a crystal whose data merge at 0.995). The refined pass is sent back only where it merges more unique reflections than its cell can hold, where it measured decisively less signal (10 %, and only where it holds no more reflections either — a wider integration disk pulls weak reflections in and dilutes the strong fraction without measuring less), or where it lost the axial rows the systematic absences are read off. One index-time veto remains and is keyed to pass 1's **lattice** rather than its group: a centred pass-1 lattice against a primitive pass-2 one. Only pass 2 is written, as the canonical `_*` output. Pass 1's merge exists to give the guard something to judge pass 2 against, so it stops short of the parts of the merge that only fill in a file — the correction surfaces, the twinning and radiation-damage analyses, the R-free flags and the amplitudes — and writes no merged files of its own. diff --git a/rugnux/Rugnux.cpp b/rugnux/Rugnux.cpp index fd6207325..68b0a7380 100644 --- a/rugnux/Rugnux.cpp +++ b/rugnux/Rugnux.cpp @@ -2505,9 +2505,10 @@ ProcessResult Rugnux::RunAllPasses(RugnuxObserver *observer) { // geometry it did not come from, a wrong header does, and a pass that lost the crystal is thrown // out by the quality guard below, which returns the run to the header geometry. // - // What starts the walk is a move of more than one step: that is the fit saying the geometry is - // somewhere else entirely rather than a fraction of a percent away, and it leaves every run - // whose fit settles beside its header at the two passes it always had. Once the run IS walking + // What starts the walk outright is a move of more than one step: that is the fit saying the + // geometry is somewhere else entirely rather than a fraction of a percent away. A smaller move + // has to be preferred by the validation frames first (below), so a run whose fit settles beside + // the geometry in hand costs two indexing probes and no pass. Once the run IS walking // it keeps walking because the steps do not get smaller in proportion to what is left: measured // on a header 3.9 % long, the first round took a third of the error and each of the next took a // third of the rest. @@ -2525,6 +2526,53 @@ ProcessResult Rugnux::RunAllPasses(RugnuxObserver *observer) { // nothing ratifies it. A walk still paying when the rounds run out keeps its best round - the // last - and says it did not converge; the quality guard below judges it against the header // pass as it judges every refined pass. + // + // The residual is not the only thing a round realises. It is the mean over every reflection, + // measured by centroids taken inside a disc centred on the prediction - so they follow the + // prediction part of the way - and dominated by the low-resolution reflections, where a distance + // and the cell scale that compensates it move every spot alike. What does see such an error is + // the re-indexing: the share of the validation spots that land on the lattice, which the + // high-resolution spots lose first. Measured on a crystal whose canonical pass ran 0.6 % long, + // the round at the fit's distance put 70.4 % of the validation spots on the lattice against + // 58.1 %, while reading the residual only 0.74 of its noise lower - and the run the residual + // stopped kept the long distance and lost 0.05 A of resolution. So a round is kept when EITHER + // realises a gain beyond its noise: the residual (HeldOutResidualFell), or the validation + // frames (ValidationEvidencePrefers). + // + // And a move of less than one step no longer settles the geometry by being small: it is asked + // of the validation frames too, but cheaply first - an indexing probe at the fit's geometry + // against one at the geometry in hand, each scoring the lattice and stopping there - and only a + // move the probes prefer pays for a re-integrated round. A move of a step or more is run as a + // round directly, as it always was, and so is every round of the walk it starts. + const auto on_lattice_pct = [](const ValidationSpotEvidence &e) { + return 100.0 * static_cast(e.on_lattice - e.by_chance) + / static_cast(std::max(1, e.spots)); + }; + // An indexing probe at a geometry: the validation evidence of the lattice indexed there, and + // nothing else - the pass stops once the lattice is scored, and the experiment is put back. + const auto index_at_geometry = [&](const std::array &g) { + const DiffractionExperiment before_probe = experiment_; + const bool probe_was = postrefine_probe_; + const bool searched_before_probe = beam_center_searched_; + beam_center_searched_ = true; // a probe scores the geometry it is given, nothing else + set_geometry(g); + indexing_probe_only_ = true; + postrefine_probe_ = false; + ProcessResult r; + try { + r = RunPipeline(observer, /*write_output=*/false, /*geometry_prepass=*/false); + } catch (const std::exception &e) { + if (IsFatalResourceError(e)) throw; + logger.Info("Two-pass: the indexing probe at distance {:.3f} mm did not complete ({})", + g[2], e.what()); + } + indexing_probe_only_ = false; + postrefine_probe_ = probe_was; + beam_center_searched_ = searched_before_probe; + experiment_ = before_probe; + ++arm_passes; + return r.validation_evidence; + }; constexpr int MAX_GEOMETRY_ROUNDS = 8; int geometry_rounds = 0; // passes run at a geometry the walk moved to int walk_passes = 0; // every pass the walk ran, the return to its best round included @@ -2537,11 +2585,14 @@ ProcessResult Rugnux::RunAllPasses(RugnuxObserver *observer) { std::string walk_rounds; // distance / realised residual per round, for the pass decision if (pass2.post_refine) { best_fit = *pass2.post_refine; - walk_rounds = fmt::format("{:.3f} mm / {:.3e}", best_geometry[2], best_fit.held_out_before); + walk_rounds = fmt::format("{:.3f} mm / {:.3e} / {:.1f}%", best_geometry[2], best_fit.held_out_before, + on_lattice_pct(pass2.validation_evidence)); } + ValidationSpotEvidence best_evidence = pass2.validation_evidence; + std::optional probe_at_best; // the indexing probe at best_geometry, once run + bool large_walk = false; // a round of this walk was started by a move of a step or more while (!cancelled_ && prepass_detector_geometry_ - && pass2.post_refine && pass2.post_refine->detector_refined - && (pass2.post_refine->large_move || geometry_rounds > 0)) { + && pass2.post_refine && pass2.post_refine->detector_refined) { if (geometry_rounds == MAX_GEOMETRY_ROUNDS) { geometry_not_converged = true; walk_stop = fmt::format("the fit still moved the geometry after {} rounds, every one of " @@ -2552,6 +2603,23 @@ ProcessResult Rugnux::RunAllPasses(RugnuxObserver *observer) { } const PostRefineResult fit = *pass2.post_refine; const std::array g = *prepass_detector_geometry_; + large_walk = large_walk || fit.large_move; + if (!large_walk) { + if (!probe_at_best) + probe_at_best = index_at_geometry(best_geometry); + const ValidationSpotEvidence probe_at_fit = index_at_geometry(g); + const bool prefers = ValidationEvidencePrefers(*probe_at_best, probe_at_fit); + logger.Info("Two-pass: the post-refinement moves the geometry by less than a step (distance " + "{:.3f} -> {:.3f} mm, beam {:.2f},{:.2f} -> {:.2f},{:.2f} px); indexing probes put " + "{:.1f}% of the validation spots on the lattice beyond chance there against {:.1f}% " + "at the geometry in hand - {}", fit.distance_before_mm, fit.distance_after_mm, + fit.beam_x_before_px, fit.beam_y_before_px, fit.beam_x_after_px, fit.beam_y_after_px, + on_lattice_pct(probe_at_fit), on_lattice_pct(*probe_at_best), + prefers ? "worth a round" : "the geometry in hand stands"); + if (!prefers) + break; + probe_at_best = probe_at_fit; + } logger.Info("Two-pass: the post-refinement at the adopted geometry still moves it " "(distance {:.3f} -> {:.3f} mm) - re-integrating and re-indexing at what it " "asks for (round {} of at most {})", fit.distance_before_mm, @@ -2561,20 +2629,26 @@ ProcessResult Rugnux::RunAllPasses(RugnuxObserver *observer) { ++geometry_rounds; ++walk_passes; const double realised = pass2.post_refine ? pass2.post_refine->held_out_before : NAN; - walk_rounds += fmt::format(", {:.3f} mm / {:.3e}", g[2], realised); - if (!pass2.post_refine || !HeldOutResidualFell(best_fit, *pass2.post_refine)) { + walk_rounds += fmt::format(", {:.3f} mm / {:.3e} / {:.1f}%", g[2], realised, + on_lattice_pct(pass2.validation_evidence)); + const bool residual_fell = pass2.post_refine && HeldOutResidualFell(best_fit, *pass2.post_refine); + const bool evidence_prefers = ValidationEvidencePrefers(best_evidence, pass2.validation_evidence); + if (!pass2.post_refine || !(residual_fell || evidence_prefers)) { walk_stop = fmt::format( - "round {} did not lower the realised held-out residual by more than its noise " - "({:.3e} against {:.3e}, standard errors {:.1e} and {:.1e})", geometry_rounds, - realised, best_fit.held_out_before, + "round {} lowered neither the realised held-out residual by more than its noise " + "({:.3e} against {:.3e}, standard errors {:.1e} and {:.1e}) nor the share of " + "validation spots off the lattice ({:.1f}% on it beyond chance against {:.1f}%)", + geometry_rounds, realised, best_fit.held_out_before, pass2.post_refine ? pass2.post_refine->held_out_before_se : NAN, - best_fit.held_out_before_se); + best_fit.held_out_before_se, on_lattice_pct(pass2.validation_evidence), + on_lattice_pct(best_evidence)); logger.Info("Two-pass: {}", walk_stop); break; } best_round = geometry_rounds; best_geometry = g; best_fit = *pass2.post_refine; + best_evidence = pass2.validation_evidence; const int pushback = ReindexPushesCellBack(fit, *pass2.post_refine); if (pushback != 0 && pushback == last_pushback) { walk_stop = fmt::format( @@ -2615,7 +2689,7 @@ ProcessResult Rugnux::RunAllPasses(RugnuxObserver *observer) { if (geometry_rounds > 0) pass2.pass_decision = fmt::format( "post-refined geometry walked {} round{} and kept round {}, detector distance {:.3f} mm " - "against the header's {:.3f}: {}; distance / realised held-out residual by round: {}", + "against the header's {:.3f}: {}; distance / realised held-out residual / validation spots on the lattice beyond chance by round: {}", geometry_rounds, geometry_rounds == 1 ? "" : "s", best_round, best_geometry[2], pass1.post_refine ? pass1.post_refine->distance_before_mm : 0.0, walk_stop.empty() ? "the fit at the kept geometry commits no further move" : walk_stop, From 6c9205a6701bfffe85c860998f2b470ffeadae94 Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Fri, 25 Sep 2026 08:57:28 +0200 Subject: [PATCH 2/4] Report the strongest index-2 superstructure class and the cell it doubles to The report-only supercell probe now reaches the results report: SUPERCELL_CLASS (parity of the primitive indices), its occupancy and Bragg-like (rocking) part against the lattice's own reflections, its , and SUPERCELL_DOUBLED_CELL, the Niggli-reduced cell to give with -C. No decision is taken on it: correct cells whose half-integer class is diffuse, and cells whose depositor kept the sub-cell, cannot be told from a real doubling by the data alone. On a crystal with a real doubled axis the reported cell matched the deposited one, and processing on it with -C brought R_free against the deposited model from 0.60 to 0.29. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01D1G8gJVAy6gp1K5Dz3NE5C --- docs/CHANGELOG.md | 1 + docs/RUGNUX_REPORT.md | 19 +++++++++++++++++++ image_analysis/IndexAndRefine.cpp | 7 +++++++ image_analysis/IndexAndRefine.h | 4 ++++ rugnux/ResultReport.cpp | 26 ++++++++++++++++++++++++++ rugnux/Rugnux.cpp | 29 +++++++++++++++++++++++++++++ rugnux/Rugnux.h | 10 ++++++++++ 7 files changed, 96 insertions(+) diff --git a/docs/CHANGELOG.md b/docs/CHANGELOG.md index 1ad3e1b3e..c88ec9e35 100644 --- a/docs/CHANGELOG.md +++ b/docs/CHANGELOG.md @@ -9,6 +9,7 @@ * Rugnux accepts a detector-modulation or absorption correction surface when held-out data support it on Fisher's z, and fits them in the order modulation, time, goniometer frame; this recovers corrections that were wrongly refused. * Rugnux settles a disagreement between the file's and the measured beam centre over an axis harmonic in either direction, so a file centre off along the spindle no longer leaves the run on a doubled axis. * Rugnux commits a post-refined geometry only when the held-out residual falls by more than its own noise, instead of by a fixed 2 %. +* Rugnux reports the strongest index-2 superstructure class of a rotation run (`SUPERCELL_CLASS` and related keys), with the doubled cell to give with `-C`; it does not change the lattice. * Rugnux also reports `R_MEAS_WEIGHTED`, the R_meas with each observation weighted as the merge weights it, so weak frames kept in the merge at low weight can be told from real disagreement; `R_MEAS` is unchanged. * Rugnux no longer stops per-frame scaling as unsettled because of a single frame that is about to be dropped as blank. * `rugnux --polarization` is documented as the polarization degree (XDS `FRACTION_OF_POLARIZATION` = (1 + p)/2); the default 0.99 suits undulators. diff --git a/docs/RUGNUX_REPORT.md b/docs/RUGNUX_REPORT.md index 4d09807a1..ad966ea99 100644 --- a/docs/RUGNUX_REPORT.md +++ b/docs/RUGNUX_REPORT.md @@ -415,6 +415,25 @@ Rings are detected far more often than they are excluded: exclusion happens only found no usable lattice. The measurement is described in [CPU/GPU data analysis ▸ Resolution and ice-ring handling](CPU_DATA_ANALYSIS_IMAGE.md#33-resolution-and-ice-ring-handling). +## Index-2 superstructure + +On rotation data, section 4 reports whether the lattice the run adopted has intensity at half-integer +positions it does not index. On 60 frames spread over the sweep, after each frame's own integration, the +lattice is predicted doubled along all three primitive axes and integrated to 3 Å; the reflections split +into eight parity classes of h, k and l, of which `0 0 0` is the lattice itself and each of the other +seven is one index-2 superstructure. Nothing is decided on it: whether a superstructure belongs in the +cell is as much the depositor's call as the data's, and several crystals whose accepted cell is the +sub-cell carry one. + +**`SUPERCELL_CLASS=`** is the parity class with the most intensity over 20–3 Å. +**`SUPERCELL_OCCUPANCY_PCT=`** is its mean intensity against the lattice's own reflections on the same +frames, and **`SUPERCELL_ROCK_PCT=`** (± **`SUPERCELL_ROCK_SE_PCT=`**) the part of it that follows the +partiality the way a Bragg reflection does, from a fit I = a + b p over the class, on the same scale. +**`SUPERCELL_I_OVER_SIGMA=`** is its mean I/σ. A class near zero on both is empty. One that is occupied +and rocks is a superstructure whose reflections this run did not integrate; **`SUPERCELL_DOUBLED_CELL=`** +is the Niggli-reduced cell the lattice would double to, to give with `-C` to process on it. One that is +occupied but hardly rocks is diffuse or disordered intensity rather than Bragg reflections. + ## Translational pseudo-symmetry Two copies of the contents of the asymmetric unit related by a pure translation that is not a lattice diff --git a/image_analysis/IndexAndRefine.cpp b/image_analysis/IndexAndRefine.cpp index f4fca96c5..8dea3f1c9 100644 --- a/image_analysis/IndexAndRefine.cpp +++ b/image_analysis/IndexAndRefine.cpp @@ -640,6 +640,7 @@ void IndexAndRefine::ProbeSupercell(std::vector image_numbers) { std::sort(image_numbers.begin(), image_numbers.end()); supercell_probe_frames_ = std::move(image_numbers); supercell_probe_ = {}; + supercell_probe_primitive_.reset(); } SupercellProbe IndexAndRefine::GetSupercellProbe() { @@ -647,6 +648,11 @@ SupercellProbe IndexAndRefine::GetSupercellProbe() { return supercell_probe_; } +std::optional IndexAndRefine::GetSupercellProbePrimitive() { + const std::unique_lock ul(supercell_probe_mutex_); + return supercell_probe_primitive_; +} + void IndexAndRefine::ProbeSupercellFrame(const DataMessage &msg, BraggPrediction &prediction, const BraggIntegrateFn &integrate, const IndexingOutcome &outcome, const BraggPredictionSettings &settings) { @@ -684,6 +690,7 @@ void IndexAndRefine::ProbeSupercellFrame(const DataMessage &msg, BraggPrediction for (int i = 0; i < 8; i++) for (int j = 0; j < 2; j++) supercell_probe_[i][j] += frame[i][j]; + supercell_probe_primitive_ = prim; } std::optional diff --git a/image_analysis/IndexAndRefine.h b/image_analysis/IndexAndRefine.h index 7d5537fe0..f3b30055c 100644 --- a/image_analysis/IndexAndRefine.h +++ b/image_analysis/IndexAndRefine.h @@ -123,6 +123,7 @@ class IndexAndRefine { // Supercell probe: the image numbers it reads, and what it has read (see ProbeSupercell). std::vector supercell_probe_frames_; SupercellProbe supercell_probe_{}; + std::optional supercell_probe_primitive_; std::mutex supercell_probe_mutex_; void ProbeSupercellFrame(const DataMessage &msg, BraggPrediction &prediction, const BraggIntegrateFn &integrate, const IndexingOutcome &outcome, @@ -189,6 +190,9 @@ public: // The probe's integrations are the caller's to keep out of anything else it counts. void ProbeSupercell(std::vector image_numbers); [[nodiscard]] SupercellProbe GetSupercellProbe(); + // The primitive lattice the parity classes are indexed in (the last probed frame's; only its metric + // is meant to be read), so a class can be turned into the cell it would double to. + [[nodiscard]] std::optional GetSupercellProbePrimitive(); // Index a single frame (no integration) with the current forced rotation lattice; used to score // first-pass sampling schemes on the real per-image path. Returns whether the frame indexed. bool IndexFrameOnly(DataMessage &msg, const SpotFindingSettings &settings); diff --git a/rugnux/ResultReport.cpp b/rugnux/ResultReport.cpp index 5f485d3de..5ed23090a 100644 --- a/rugnux/ResultReport.cpp +++ b/rugnux/ResultReport.cpp @@ -1075,6 +1075,32 @@ ReportDocument BuildReportDocument(const std::string &output_prefix, } } + // ---- index-2 superstructure + // Report only, like the harmonic above: on 60 frames the probe integrates the lattice doubled + // along all three primitive axes and reads the seven half-integer parity classes. Whether a + // class that is there belongs in the cell is the depositor's call as much as the data's - + // several crystals whose accepted cell is the sub-cell carry one - so no verdict is printed. + if (result.supercell) { + const auto &sc = *result.supercell; + Add(s, KeyText("SUPERCELL_CLASS", fmt::format("{} {} {}", sc.h, sc.k, sc.l))); + Add(s, KeyReal("SUPERCELL_OCCUPANCY_PCT", sc.occupancy_pct, "{:.1f}")); + Add(s, KeyReal("SUPERCELL_ROCK_PCT", sc.rock_pct, "{:.1f}")); + Add(s, KeyReal("SUPERCELL_ROCK_SE_PCT", sc.rock_se_pct, "{:.1f}")); + Add(s, KeyReal("SUPERCELL_I_OVER_SIGMA", sc.i_over_sigma, "{:.2f}")); + Add(s, KeyText("SUPERCELL_DOUBLED_CELL", CellString(sc.doubled_cell))); + Add(s, Blank()); + Add(s, Prose( + " SUPERCELL_CLASS is the half-integer class (parity of h, k, l in the primitive cell of the\n" + " lattice) with the most intensity over 20-3 A, read on 60 frames. OCCUPANCY is its mean\n" + " intensity against the lattice's own reflections; ROCK is the part of it that follows the\n" + " partiality the way a Bragg reflection does, on the same scale. A class near zero on both is\n" + " empty. One that is occupied and rocks is an index-2 superstructure the lattice does not\n" + " index: its reflections are not integrated. To process on the doubled lattice, give\n" + " SUPERCELL_DOUBLED_CELL with -C. An occupied class that does not rock is diffuse or\n" + " disordered intensity rather than Bragg reflections.")); + Add(s, Blank()); + } + // ---- powder contamination // A crystalline phase other than the crystal, diffracting as rings among its reflections. // Measured in the pre-scan on every run, from the spots found there for the width and the diff --git a/rugnux/Rugnux.cpp b/rugnux/Rugnux.cpp index fd6207325..55b6918f4 100644 --- a/rugnux/Rugnux.cpp +++ b/rugnux/Rugnux.cpp @@ -5312,6 +5312,35 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b } logger.Info("Supercell probe (the eight parity classes of the doubled primitive cell; , " "mean intensity, and the fit I = a + b p against the lattice's own):{}", table); + + // The class with the most intensity over 20-3 A goes into the report, with the cell it would + // double the lattice to: the vectors t of the primitive lattice with (h,k,l).t even - twice + // one axis the class is odd on, and each other axis plus that one where the class is odd + // on it too - Niggli-reduced. + const auto band = [&](int c) { SupercellProbeClass x = probe[c][0]; x += probe[c][1]; return x; }; + const SupercellProbeFit f0 = FitSupercellProbe(band(0)); + const auto prim = indexer->GetSupercellProbePrimitive(); + int best = 0; + for (int c = 1; c < 8; c++) + if (best == 0 || FitSupercellProbe(band(c)).mean_i > FitSupercellProbe(band(best)).mean_i) + best = c; + if (prim && f0.mean_i > 0.0 && f0.b != 0.0) { + const SupercellProbeFit f = FitSupercellProbe(band(best)); + const int p[3] = {(best >> 2) & 1, (best >> 1) & 1, best & 1}; + const Coord v[3] = {prim->Vec0(), prim->Vec1(), prim->Vec2()}; + const int i = p[0] ? 0 : p[1] ? 1 : 2; + Coord t[3]; + for (int j = 0; j < 3; j++) + t[j] = j == i ? v[i] * 2.0f : v[j] + v[i] * static_cast(p[j]); + ProcessResult::SupercellClass sc; + sc.h = p[0]; sc.k = p[1]; sc.l = p[2]; + sc.occupancy_pct = 100.0 * f.mean_i / f0.mean_i; + sc.rock_pct = 100.0 * f.b / f0.b; + sc.rock_se_pct = 100.0 * f.b_se / f0.b; + sc.i_over_sigma = f.mean_i_over_sigma; + sc.doubled_cell = CrystalLattice(t[0], t[1], t[2]).NiggliReduce().GetUnitCell(); + result.supercell = sc; + } } } diff --git a/rugnux/Rugnux.h b/rugnux/Rugnux.h index 29e1c315a..77dd4c715 100644 --- a/rugnux/Rugnux.h +++ b/rugnux/Rugnux.h @@ -362,6 +362,16 @@ struct ProcessResult { // Higher-order contamination of the beam, read off the spots the lattice did not take. Measured // on the images themselves, so it is there whether or not anything merged. HarmonicContaminationResult harmonic; + // The strongest index-2 superstructure class the supercell probe read (20-3 A, 60 frames): which + // parity class of the primitive cell, its mean intensity and the part that rocks like a Bragg + // reflection, both against the lattice's own, and the cell the lattice would double to were it + // real. Report only - nothing decides on it. + struct SupercellClass { + int h = 0, k = 0, l = 0; // parity of the primitive indices + double occupancy_pct = 0.0, rock_pct = 0.0, rock_se_pct = 0.0, i_over_sigma = 0.0; + UnitCell doubled_cell{}; + }; + std::optional supercell; // Powder contamination measured in the pre-scan: a crystalline phase other than the crystal, // diffracting as rings among its reflections. Measured and reported on every run. PowderRings powder; From dd18d8c21dafb5d2ee3038ec0b3eaee96b2db008 Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Fri, 25 Sep 2026 09:36:25 +0200 Subject: [PATCH 3/4] Warn SUPERCELL_POSSIBLE where the half-integer class looks like Bragg reflections Advisory only: the class is measured ( >= 0.5) and its rocking part is at least 2% of the lattice's, three standard errors clear. The warning asks the user to process both settings - the run's cell and the reported doubled cell with -C - and compare them in refinement. Of 105 battery rotation sets it names 8: both crystals whose accepted cell is the doubled one, one pseudo-translation whose doubled description is an accepted alternative, and five correct sub-cells with weak ordered half-integer intensity. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01D1G8gJVAy6gp1K5Dz3NE5C --- docs/CHANGELOG.md | 2 +- docs/RUGNUX_REPORT.md | 9 ++++++++- rugnux/ReportDocument.h | 1 + rugnux/ResultReport.cpp | 16 ++++++++++++++++ 4 files changed, 26 insertions(+), 2 deletions(-) diff --git a/docs/CHANGELOG.md b/docs/CHANGELOG.md index c88ec9e35..9e1b51dd0 100644 --- a/docs/CHANGELOG.md +++ b/docs/CHANGELOG.md @@ -9,7 +9,7 @@ * Rugnux accepts a detector-modulation or absorption correction surface when held-out data support it on Fisher's z, and fits them in the order modulation, time, goniometer frame; this recovers corrections that were wrongly refused. * Rugnux settles a disagreement between the file's and the measured beam centre over an axis harmonic in either direction, so a file centre off along the spindle no longer leaves the run on a doubled axis. * Rugnux commits a post-refined geometry only when the held-out residual falls by more than its own noise, instead of by a fixed 2 %. -* Rugnux reports the strongest index-2 superstructure class of a rotation run (`SUPERCELL_CLASS` and related keys), with the doubled cell to give with `-C`; it does not change the lattice. +* Rugnux reports the strongest index-2 superstructure class of a rotation run (`SUPERCELL_CLASS` and related keys) with the doubled cell to give with `-C`, and warns (`SUPERCELL_POSSIBLE`) where that class looks like Bragg reflections, so both settings can be tried; it does not change the lattice. * Rugnux also reports `R_MEAS_WEIGHTED`, the R_meas with each observation weighted as the merge weights it, so weak frames kept in the merge at low weight can be told from real disagreement; `R_MEAS` is unchanged. * Rugnux no longer stops per-frame scaling as unsettled because of a single frame that is about to be dropped as blank. * `rugnux --polarization` is documented as the polarization degree (XDS `FRACTION_OF_POLARIZATION` = (1 + p)/2); the default 0.99 suits undulators. diff --git a/docs/RUGNUX_REPORT.md b/docs/RUGNUX_REPORT.md index ad966ea99..82f74bf2a 100644 --- a/docs/RUGNUX_REPORT.md +++ b/docs/RUGNUX_REPORT.md @@ -148,7 +148,7 @@ and a run that produced a textbook data set read identically for their first thr `SYMMETRY_AMBIGUITY`, `CENTERING_UNTESTED`, `UNUSABLE_MERGE`, `LOW_COMPLETENESS`, `SWEEP_GAPS`, `GONIO_SCALE`, `SPINDLE_CAP`, `ANISOTROPY`, `TWINNING`, `PSEUDO_TRANSLATION`, `LATTICE_TRANSLATION`, `MODEL_HAND`, `MODEL_NOT_VALIDATED`, `CANCELLED`, `RESOLUTION_FIT`, - `FLIGHT_PATH`, `GEOMETRY_NOT_CONVERGED`. `NONE` when nothing fired. A code appears if and only if its + `FLIGHT_PATH`, `GEOMETRY_NOT_CONVERGED`, `SUPERCELL_POSSIBLE`. `NONE` when nothing fired. A code appears if and only if its warning fired, so the flags and the `WARNING:` lines are two renderings of one list — the closed type for machinery, the open sentence for a person. - **`WARNING_COUNT=`** and the **`WARNING:`** lines follow, in the same section. They are what they @@ -434,6 +434,13 @@ and rocks is a superstructure whose reflections this run did not integrate; **`S is the Niggli-reduced cell the lattice would double to, to give with `-C` to process on it. One that is occupied but hardly rocks is diffuse or disordered intensity rather than Bragg reflections. +**`SUPERCELL_POSSIBLE=`** is `TRUE`, with a warning under the `SUPERCELL_POSSIBLE` flag, where the +class is measured (`SUPERCELL_I_OVER_SIGMA=` at least 0.5) and part of it rocks like Bragg reflections +(`SUPERCELL_ROCK_PCT=` at least 2 %, three standard errors clear of zero). It is advice to check, not a +finding: the same numbers come from a real doubled cell and from a correct cell with weak ordered +intensity between its reflections, and which of the two a structure is decided by refinement. Process +both settings - the run's cell, and `SUPERCELL_DOUBLED_CELL=` given with `-C` - and compare them there. + ## Translational pseudo-symmetry Two copies of the contents of the asymmetric unit related by a pure translation that is not a lattice diff --git a/rugnux/ReportDocument.h b/rugnux/ReportDocument.h index be8690150..3b72acc0f 100644 --- a/rugnux/ReportDocument.h +++ b/rugnux/ReportDocument.h @@ -43,6 +43,7 @@ namespace PathologyCode { // in how exact it is. ANISOTROPY and RESOLUTION_FIT are already shared this way. constexpr const char *LATTICE_TRANSLATION = "LATTICE_TRANSLATION"; // a translation the cell does not declare constexpr const char *HARMONIC_CONTAMINATION = "HARMONIC_CONTAMINATION"; // the beam carries a higher harmonic + constexpr const char *SUPERCELL_POSSIBLE = "SUPERCELL_POSSIBLE"; // Bragg-like intensity at half-integer positions constexpr const char *GEOMETRY_NOT_CONVERGED = "GEOMETRY_NOT_CONVERGED"; // the geometry walk ran out of rounds constexpr const char *SCALING_NOT_CONVERGED = "SCALING_NOT_CONVERGED"; // the per-frame scales hit their cap still moving } diff --git a/rugnux/ResultReport.cpp b/rugnux/ResultReport.cpp index 5ed23090a..2b64d39ad 100644 --- a/rugnux/ResultReport.cpp +++ b/rugnux/ResultReport.cpp @@ -1088,6 +1088,22 @@ ReportDocument BuildReportDocument(const std::string &output_prefix, Add(s, KeyReal("SUPERCELL_ROCK_SE_PCT", sc.rock_se_pct, "{:.1f}")); Add(s, KeyReal("SUPERCELL_I_OVER_SIGMA", sc.i_over_sigma, "{:.2f}")); Add(s, KeyText("SUPERCELL_DOUBLED_CELL", CellString(sc.doubled_cell))); + // Advisory, not a decision: the class is measured ( at least 0.5) and part of it + // rocks like Bragg reflections (at least 2% of the lattice's, three standard errors clear). + // On the battery's 105 rotation sets this names 8 - both crystals whose accepted cell is the + // doubled one, a pseudo-translation whose doubled description is an accepted alternative, + // and five correct sub-cells with real but weak half-integer intensity. + const bool possible = sc.i_over_sigma >= 0.5 && sc.rock_pct >= 2.0 && sc.rock_pct > 3.0 * sc.rock_se_pct; + Add(s, KeyBool("SUPERCELL_POSSIBLE", possible)); + if (possible) + Warn(doc, PathologyCode::SUPERCELL_POSSIBLE, + fmt::format("Half-integer reflections of parity {} {} {} carry Bragg-like intensity ({:.0f}% " + "of the lattice's, rocking part {:.1f}%): the true cell may be twice as large. " + "Check carefully, preferably by processing both settings - this cell and " + "-C \"{:.2f},{:.2f},{:.2f},{:.2f},{:.2f},{:.2f}\" - and comparing them in refinement", + sc.h, sc.k, sc.l, sc.occupancy_pct, sc.rock_pct, + sc.doubled_cell.a, sc.doubled_cell.b, sc.doubled_cell.c, + sc.doubled_cell.alpha, sc.doubled_cell.beta, sc.doubled_cell.gamma)); Add(s, Blank()); Add(s, Prose( " SUPERCELL_CLASS is the half-integer class (parity of h, k, l in the primitive cell of the\n" From 3de894989c94a0cd6314acab3b0ee81e272e10ef Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Fri, 25 Sep 2026 10:20:28 +0200 Subject: [PATCH 4/4] Report summary: a Supercell line beside pseudo-symmetry and twinning The SUPERCELL_POSSIBLE advice now also sits in the summary at the top of the report, with the doubled cell ready for -C; the rule is shared with the warning. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01D1G8gJVAy6gp1K5Dz3NE5C --- docs/RUGNUX_REPORT.md | 4 +++- rugnux/ResultReport.cpp | 26 ++++++++++++++++++++------ 2 files changed, 23 insertions(+), 7 deletions(-) diff --git a/docs/RUGNUX_REPORT.md b/docs/RUGNUX_REPORT.md index 82f74bf2a..2c47cdf8c 100644 --- a/docs/RUGNUX_REPORT.md +++ b/docs/RUGNUX_REPORT.md @@ -439,7 +439,9 @@ class is measured (`SUPERCELL_I_OVER_SIGMA=` at least 0.5) and part of it rocks (`SUPERCELL_ROCK_PCT=` at least 2 %, three standard errors clear of zero). It is advice to check, not a finding: the same numbers come from a real doubled cell and from a correct cell with weak ordered intensity between its reflections, and which of the two a structure is decided by refinement. Process -both settings - the run's cell, and `SUPERCELL_DOUBLED_CELL=` given with `-C` - and compare them there. +both settings - the run's cell, and `SUPERCELL_DOUBLED_CELL=` given with `-C` - and compare them there. The +summary at the top of the report carries the same advice on its `Supercell` line, beside +pseudo-symmetry and twinning. ## Translational pseudo-symmetry diff --git a/rugnux/ResultReport.cpp b/rugnux/ResultReport.cpp index 2b64d39ad..f56b5e7cb 100644 --- a/rugnux/ResultReport.cpp +++ b/rugnux/ResultReport.cpp @@ -50,6 +50,15 @@ namespace { const char *BANNER = " ******************************************************************************"; + // Advisory, not a decision: the half-integer class is measured ( at least 0.5) and part of + // it rocks like Bragg reflections (at least 2% of the lattice's, three standard errors clear). On + // the battery's 105 rotation sets this names 8 - both crystals whose accepted cell is the doubled + // one, a pseudo-translation whose doubled description is an accepted alternative, and five correct + // sub-cells with real but weak half-integer intensity. + bool SupercellPossible(const ProcessResult::SupercellClass &sc) { + return sc.i_over_sigma >= 0.5 && sc.rock_pct >= 2.0 && sc.rock_pct > 3.0 * sc.rock_se_pct; + } + std::string CellString(const UnitCell &c) { return fmt::format("{:.3f} {:.3f} {:.3f} {:.3f} {:.3f} {:.3f}", c.a, c.b, c.c, c.alpha, c.beta, c.gamma); @@ -1088,12 +1097,7 @@ ReportDocument BuildReportDocument(const std::string &output_prefix, Add(s, KeyReal("SUPERCELL_ROCK_SE_PCT", sc.rock_se_pct, "{:.1f}")); Add(s, KeyReal("SUPERCELL_I_OVER_SIGMA", sc.i_over_sigma, "{:.2f}")); Add(s, KeyText("SUPERCELL_DOUBLED_CELL", CellString(sc.doubled_cell))); - // Advisory, not a decision: the class is measured ( at least 0.5) and part of it - // rocks like Bragg reflections (at least 2% of the lattice's, three standard errors clear). - // On the battery's 105 rotation sets this names 8 - both crystals whose accepted cell is the - // doubled one, a pseudo-translation whose doubled description is an accepted alternative, - // and five correct sub-cells with real but weak half-integer intensity. - const bool possible = sc.i_over_sigma >= 0.5 && sc.rock_pct >= 2.0 && sc.rock_pct > 3.0 * sc.rock_se_pct; + const bool possible = SupercellPossible(sc); Add(s, KeyBool("SUPERCELL_POSSIBLE", possible)); if (possible) Warn(doc, PathologyCode::SUPERCELL_POSSIBLE, @@ -2051,6 +2055,16 @@ ReportDocument BuildReportDocument(const std::string &output_prefix, : std::string("no indication")); if (result.twinning.l_test_pairs > 0) row("Twinning", TwinningVerdictLine(result.twinning)); + if (result.supercell) + row("Supercell", SupercellPossible(*result.supercell) + ? fmt::format("POSSIBLE - class {} {} {} rocks like Bragg reflections ({:.1f}%); " + "process also with -C \"{:.2f},{:.2f},{:.2f},{:.2f},{:.2f},{:.2f}\" " + "and compare", result.supercell->h, result.supercell->k, + result.supercell->l, result.supercell->rock_pct, + result.supercell->doubled_cell.a, result.supercell->doubled_cell.b, + result.supercell->doubled_cell.c, result.supercell->doubled_cell.alpha, + result.supercell->doubled_cell.beta, result.supercell->doubled_cell.gamma) + : std::string("no indication")); const double db = result.merge_statistics.radiation_damage_delta_b; if (std::isfinite(db)) row("Radiation damage", fmt::format("relative B {:+.2f} A^2 over the sweep", db));