From 0ca159449f71f7767f0b90c88d9627381ca50443 Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Fri, 31 Jul 2026 14:03:34 +0200 Subject: [PATCH] Bragg integration: integrate as far as the detector reaches, not to a fixed 1.0 A BraggIntegrationSettings::DMinLimit_A had a setter that nothing anywhere called, so it was always its 1.0 A default - in rugnux, the viewer and the broker alike, with no option or API field to change it. It feeds the predictor as high_res_A, which discards any reflection with |q| > 1/d_min, so integration simply stopped at 1.0 A however far the detector reached. Five of the 33 rotation test datasets have detectors reaching past it, down to 0.981 A. On one of them, run with no resolution limit, the shell table ended dead at 1.00 A with that shell still at CC1/2 55.6% and 3.4 - cut mid-shell rather than fading out. This branch had already made the sibling limits detector-driven (spot finding, scaling), so the pipeline was finding spots the detector could see and then refusing to integrate them. Make it a std::optional: unset means as far as the detector reaches, a value limits. The limit is only a bound on how far the lattice walk goes, never a second opinion on what is measurable - both predictors independently drop reflections that miss the detector (BraggPrediction.cpp, BraggPredictionRot.cpp) - which is what makes the detector's own reach the right default. rugnux gains --integration-high-resolution (0 = no limit, as for --spot-high-resolution); the derived per-axis prediction range resolves against the same number, so the two cannot drift. Full battery: 30/33 space groups, unchanged from before, 0 failures and the same three known mismatches; 22 of 32 crystals bit-identical and nothing worse than 5 observations in ~500k. The datasets that gain do so because their detector reached past 1.0 A - the effect is understated here because the harness caps each merge at the XDS resolution anyway. Co-Authored-By: Claude Opus 5 (1M context) --- common/BraggIntegrationSettings.cpp | 13 +++++++------ common/BraggIntegrationSettings.h | 11 ++++++++--- docs/CHANGELOG.md | 1 + docs/RUGNUX.md | 1 + image_analysis/IndexAndRefine.cpp | 25 ++++++++++++++++++------- rugnux/RugnuxCommandLine.cpp | 4 ++++ rugnux/rugnux_cli.cpp | 18 +++++++++++++++++- 7 files changed, 56 insertions(+), 17 deletions(-) diff --git a/common/BraggIntegrationSettings.cpp b/common/BraggIntegrationSettings.cpp index fc2c8aae..00f0fb0d 100644 --- a/common/BraggIntegrationSettings.cpp +++ b/common/BraggIntegrationSettings.cpp @@ -44,11 +44,12 @@ BraggIntegrationSettings &BraggIntegrationSettings::R3(float input) { return *this; } -BraggIntegrationSettings &BraggIntegrationSettings::DMinLimit_A(float input) { - check_finite("Minimum d-spacing", input); - check_min("Minimum d-spacing", input, 0.5); - check_max("Minimum d-spacing", input, 100.0); - +BraggIntegrationSettings &BraggIntegrationSettings::DMinLimit_A(std::optional input) { + if (input) { + check_finite("Minimum d-spacing", *input); + check_min("Minimum d-spacing", *input, 0.5); + check_max("Minimum d-spacing", *input, 100.0); + } d_min_limit_A = input; return *this; } @@ -89,7 +90,7 @@ float BraggIntegrationSettings::GetR3() const { return r_3; } -float BraggIntegrationSettings::GetDMinLimit_A() const { +std::optional BraggIntegrationSettings::GetDMinLimit_A() const { return d_min_limit_A; } diff --git a/common/BraggIntegrationSettings.h b/common/BraggIntegrationSettings.h index 22ce36b6..2e3e5b23 100644 --- a/common/BraggIntegrationSettings.h +++ b/common/BraggIntegrationSettings.h @@ -17,7 +17,12 @@ class BraggIntegrationSettings { float r_1 = 4; float r_2 = 6; float r_3 = 10; - float d_min_limit_A = 1.0; + // Integration/prediction resolution limit. Unset means "as far as the detector reaches", resolved + // from the geometry where it is used. The predictor independently rejects any reflection that misses + // the detector, so this is a bound on how far the lattice walk goes rather than a second opinion on + // what is measurable - a fixed default simply truncated every experiment whose detector reached + // past it. + std::optional d_min_limit_A; std::optional fixed_profile_radius; float minimum_sigma_in_regards_to_i = 0.02; // Symmetric trimmed-mean fraction for the r2..r3 background ring: drop the lowest and highest this @@ -40,7 +45,7 @@ public: BraggIntegrationSettings& R1(float input); BraggIntegrationSettings& R2(float input); BraggIntegrationSettings& R3(float input); - BraggIntegrationSettings& DMinLimit_A(float input); + BraggIntegrationSettings& DMinLimit_A(std::optional input); BraggIntegrationSettings& FixedProfileRadius_recipA(std::optional input); BraggIntegrationSettings& Integrator(IntegratorMode input); BraggIntegrationSettings& BackgroundTrimFraction(float input); @@ -52,7 +57,7 @@ public: [[nodiscard]] float GetR2() const; [[nodiscard]] float GetR3() const; [[nodiscard]] std::optional GetFixedProfileRadius_recipA() const; - [[nodiscard]] float GetDMinLimit_A() const; + [[nodiscard]] std::optional GetDMinLimit_A() const; [[nodiscard]] float GetMinimumSigmaInRegardsToI() const; [[nodiscard]] float GetBackgroundTrimFraction() const; diff --git a/docs/CHANGELOG.md b/docs/CHANGELOG.md index 2d4de63f..2f1e1513 100644 --- a/docs/CHANGELOG.md +++ b/docs/CHANGELOG.md @@ -13,6 +13,7 @@ This is an UNSTABLE release. It includes many experimental features, as well as * rugnux: Per-image **geometry refinement** weights spots by confidence. * rugnux: `--min-image-cc` now works for rotation data (opt-in). * rugnux: New `--search-min-zeta` (rotation default 0.85). +* Bragg integration: the prediction/integration resolution limit now defaults to **as far as the detector reaches**, instead of a hardcoded 1.0 Å that nothing could change. Any experiment whose detector reached past 1.0 Å was silently losing everything beyond it (5 of 33 rotation test datasets do, down to 0.981 Å); the reflections were found by spot finding and then never integrated. rugnux `--integration-high-resolution ` sets a limit by hand (0 = no limit). * Bragg prediction: how far the predictor walks the lattice is now a setting (`bragg_integration_settings.max_hkl`) rather than a hardcoded 100, for both rotation and stills. Offline (rugnux, viewer) it is **derived per crystal** from the refined cell as `ceil(max(a,b,c)/d_min) + 1` — an exact bound, so it recovers the outermost reflections that a fixed 100 was quietly truncating on a long-axis or high-resolution dataset, while costing less than before on a small cell. Measured: up to +17% observations and high-shell CC1/2 15.1% → 25.8%, with small cells bit-identical. `rugnux --max-hkl ` overrides it. Online the broker bootstraps 100, so a live acquisition keeps a predictable per-image cost whatever crystal is mounted; it is settable via the API and the frontend. * rugnux: Reports how close a symmetry axis lies to the spindle. * rugnux: Azimuthal-integration and spot-finding resolution limits default to the detector — including **rotation** data, which no longer keeps a 1.5 Å spot-finding limit (`--spot-high-resolution ` still sets one). diff --git a/docs/RUGNUX.md b/docs/RUGNUX.md index 9f46d831..11307256 100644 --- a/docs/RUGNUX.md +++ b/docs/RUGNUX.md @@ -268,6 +268,7 @@ Integration: | `--integrator ` | Spot integrator: `gaussian` (profile-fit, default) \| `empirical` \| `boxsum` (classical fallback) | | `--integration-radius ` | Signal-box radius `r1`, or `r1,r2,r3` (px). One value ⇒ `r2=r1+2`, `r3=r1+4` | | `--background-trim ` | Monochromatic (rotation + still): symmetric trimmed-mean fraction for the background ring, 0≤f<0.5 (default 0.10; 0 = plain mean) — removes the high-side bias that over-subtracts weak high-angle spots | +| `--integration-high-resolution ` | High-resolution limit for prediction and integration. Omitted (or 0) means integration extends as far as the detector reaches — which is what the predictor can place on the detector anyway, since it rejects reflections that miss it. Set a value to integrate less than the detector offers | | `--max-hkl ` | Predict reflections with \|h\|,\|k\|,\|l\| ≤ `n` (max 511). By default this is derived per crystal from the refined cell as `ceil(max(a,b,c)/d_min) + 1`, which is the exact bound: the predictor keeps only \|q\| ≤ 1/d_min and `h = a·q`, so no reflection can lie outside it and no candidate inside it is wasted on a shorter axis. Set it only to override that | | `--bandwidth ` | Relative X-ray bandwidth FWHM (e.g. `0.01` for a 1% DMM); default from file or 0 (monochromatic) | diff --git a/image_analysis/IndexAndRefine.cpp b/image_analysis/IndexAndRefine.cpp index a6655a4f..d0dbf5cb 100644 --- a/image_analysis/IndexAndRefine.cpp +++ b/image_analysis/IndexAndRefine.cpp @@ -32,22 +32,33 @@ namespace { return static_cast(std::ceil(axis_A / d_min_A)) + 1; } + // The resolution the prediction walks out to: an explicit setting, else as far as the detector + // reaches. The predictor drops any reflection that misses the detector anyway, so this bound cannot + // add reflections the geometry does not offer - it only decides how much of the lattice is examined, + // which is why the detector's own reach is the right default and a fixed number was not. A geometry + // with no scattering angle at all (no distance, no wavelength) reports 0, which is not a resolution + // to divide by; nothing can be predicted from it either way. + float IntegrationDMin_A(const DiffractionExperiment &experiment) { + if (const auto fixed = experiment.GetBraggIntegrationSettings().GetDMinLimit_A()) + return *fixed; + const float detector_A = experiment.GetDetectorMaxResolution_A(); + return detector_A > 0.0f ? detector_A : 1.0f; + } + // An explicit setting is enforced as given, on every index - it is one number, deliberately, because // it exists to bound the work rather than to describe the crystal. Otherwise the cell decides. The // online path always carries a value (the broker bootstraps one and the API can change it), so a // live acquisition never has its per-frame cost decided by whichever crystal was mounted. void ApplyPredictionRange(BraggPredictionSettings &settings, const DiffractionExperiment &experiment, const CrystalLattice &latt) { - const auto &bragg = experiment.GetBraggIntegrationSettings(); - if (const auto fixed = bragg.GetMaxHKL()) { + if (const auto fixed = experiment.GetBraggIntegrationSettings().GetMaxHKL()) { settings.max_h = settings.max_k = settings.max_l = *fixed; return; } const UnitCell cell = latt.GetUnitCell(); - const float d_min_A = bragg.GetDMinLimit_A(); - settings.max_h = MaxIndexForAxis(cell.a, d_min_A); - settings.max_k = MaxIndexForAxis(cell.b, d_min_A); - settings.max_l = MaxIndexForAxis(cell.c, d_min_A); + settings.max_h = MaxIndexForAxis(cell.a, settings.high_res_A); + settings.max_k = MaxIndexForAxis(cell.b, settings.high_res_A); + settings.max_l = MaxIndexForAxis(cell.c, settings.high_res_A); } } @@ -466,7 +477,7 @@ void IndexAndRefine::QuickPredictAndIntegrate(DataMessage &msg, }; BraggPredictionSettings settings_prediction{ - .high_res_A = experiment.GetBraggIntegrationSettings().GetDMinLimit_A(), + .high_res_A = IntegrationDMin_A(experiment), .ewald_dist_cutoff = ewald_dist_cutoff, // Centering is a hypothesis to confirm, not assume: with no user-fixed space group, predict // in P so the centering-absent reflections are integrated and the space-group search can diff --git a/rugnux/RugnuxCommandLine.cpp b/rugnux/RugnuxCommandLine.cpp index 04faa0cd..0e73f3bd 100644 --- a/rugnux/RugnuxCommandLine.cpp +++ b/rugnux/RugnuxCommandLine.cpp @@ -126,6 +126,10 @@ std::string RugnuxCommandLine(const ProcessConfig &config, add("--background-trim", num(bragg.GetBackgroundTrimFraction())); if (const auto max_hkl = bragg.GetMaxHKL()) add("--max-hkl", std::to_string(*max_hkl)); + // Unset means "to the detector edge"; emitting the resolved number would pin it to this run's + // geometry, so leave the flag out and let it resolve again. + if (const auto d_min = bragg.GetDMinLimit_A()) + add("--integration-high-resolution", num(*d_min)); if (config.rotation_indexing) { if (config.two_pass_rotation) diff --git a/rugnux/rugnux_cli.cpp b/rugnux/rugnux_cli.cpp index 65694c94..152ae0d5 100644 --- a/rugnux/rugnux_cli.cpp +++ b/rugnux/rugnux_cli.cpp @@ -117,6 +117,7 @@ void print_usage() { std::cout << " Integration" << std::endl; std::cout << " --bandwidth Relative X-ray bandwidth FWHM (e.g. 0.01 for 1% DMM); default from file or 0" << std::endl; std::cout << " --integration-radius Signal-box radius r1, or r1,r2,r3 (px). One value => r2=r1+2, r3=r1+4" << std::endl; + std::cout << " --integration-high-resolution High resolution limit for prediction/integration. If omitted (or 0), integration extends as far as the detector reaches" << std::endl; std::cout << " --max-hkl Predict reflections with |h|,|k|,|l| <= n. Default: derived per crystal from the refined cell (ceil(longest axis / d_min) + 1), which is the exact bound - set it only to override that" << std::endl; std::cout << " --background-trim Monochromatic (rotation + still): symmetric trimmed-mean fraction for the background ring (0<=f<0.5, default 0.10; 0 = plain mean). Removes the high-side bias that over-subtracts weak high-angle spots (broadband data keep the sigma-clip instead)" << std::endl; std::cout << " --integrator Spot integrator boxsum|gaussian|empirical (default: gaussian profile-fit; boxsum is the classical fallback)" << std::endl; @@ -166,6 +167,7 @@ enum { OPT_INTEGRATION_RADIUS, OPT_BACKGROUND_TRIM, OPT_MAX_HKL, + OPT_INTEGRATION_HIGH_RES, OPT_REJECT_OUTLIERS, OPT_REFERENCE_COLUMN, OPT_MODEL, @@ -274,6 +276,7 @@ static option long_options[] = { {"integration-radius", required_argument, nullptr, OPT_INTEGRATION_RADIUS}, {"background-trim", required_argument, nullptr, OPT_BACKGROUND_TRIM}, {"max-hkl", required_argument, nullptr, OPT_MAX_HKL}, + {"integration-high-resolution", required_argument, nullptr, OPT_INTEGRATION_HIGH_RES}, {"integrator", required_argument, nullptr, OPT_INTEGRATOR}, {"simple-stills", no_argument, nullptr, OPT_SIMPLE_STILLS}, {"detect-ice-rings", optional_argument, nullptr, OPT_DETECT_ICE_RINGS}, @@ -553,7 +556,8 @@ static int RunRugnux(int argc, char **argv) { std::optional report_shell_count; // --resolution-shells std::optional integration_radius_arg; // "r1" or "r1,r2,r3" std::optional background_trim_arg; // --background-trim: background-ring trimmed-mean fraction - std::optional max_hkl_arg; // --max-hkl: half-width of the predicted hkl cube + std::optional max_hkl_arg; // --max-hkl: half-width of the predicted hkl box + std::optional integration_d_min_arg; // --integration-high-resolution; unset = detector reach std::optional integrator_mode; // --integrator boxsum|gaussian|empirical bool simple_stills_flag = false; // --simple-stills: disable the default stills partiality post-refinement std::optional outlier_reject_nsigma; // merge per-observation outlier rejection @@ -848,6 +852,9 @@ static int RunRugnux(int argc, char **argv) { case OPT_MAX_HKL: max_hkl_arg = parse_number_arg(optarg, "--max-hkl", logger, 1, 511); break; + case OPT_INTEGRATION_HIGH_RES: + integration_d_min_arg = parse_double_arg(optarg, "--integration-high-resolution", logger); + break; case OPT_INTEGRATOR: if (strcmp(optarg, "boxsum") == 0) integrator_mode = IntegratorMode::BoxSum; else if (strcmp(optarg, "gaussian") == 0) integrator_mode = IntegratorMode::ProfileGaussian; @@ -1497,6 +1504,15 @@ static int RunRugnux(int argc, char **argv) { : "profile (empirical)"); } + if (integration_d_min_arg) { + BraggIntegrationSettings bis = experiment.GetBraggIntegrationSettings(); + // 0 spells "no limit" for the sibling resolution options, so it has to mean the same here. + bis.DMinLimit_A(*integration_d_min_arg > 0.0 + ? std::optional(static_cast(*integration_d_min_arg)) + : std::nullopt); + experiment.ImportBraggIntegrationSettings(bis); + } + if (max_hkl_arg) { BraggIntegrationSettings bis = experiment.GetBraggIntegrationSettings(); bis.MaxHKL(static_cast(*max_hkl_arg));