From 76e88b0fca6a9bef12414ee1a80028dddca8e73a Mon Sep 17 00:00:00 2001 From: leonarski_f Date: Sun, 28 Jun 2026 11:50:31 +0200 Subject: [PATCH] integration: sigma-clip the background ring to de-bias high-resolution stills Weak high-resolution still reflections were systematically over-subtracted: a bandwidth-streaked high-res spot (or a neighbour) leaks into the r2-r3 background annulus and biases its mean high, so the subtracted background is too large and the merged high-resolution intensities go negative (seen as reproducibly negative at 100% completeness and high multiplicity past ~1.9 A). Add one high-outlier sigma-clip pass to the box-sum background (reject ring pixels above mean + 3*sqrt(mean), recompute) so the contamination no longer inflates it. A clean Poisson background is essentially unchanged (~0.1% exceed the cut). On the HEWL serial-stills jet this de-biases the high-res band - 2.03 A 0.9 -> 1.6, 1.79 A -0.1 -> +0.7 - extends the usable resolution ~2.2 -> ~2.0 A and improves overall R-meas 130 -> 124%, with CC1/2 and CC-vs-reference neutral. The rotation crystal is unchanged (ISa 19.1), its clean backgrounds being barely clipped. Co-Authored-By: Claude Opus 4.8 --- .../bragg_integration/ProfileIntegrate2D.cpp | 14 ++++++++++++++ 1 file changed, 14 insertions(+) diff --git a/image_analysis/bragg_integration/ProfileIntegrate2D.cpp b/image_analysis/bragg_integration/ProfileIntegrate2D.cpp index c43a5d2d..95cfeed1 100644 --- a/image_analysis/bragg_integration/ProfileIntegrate2D.cpp +++ b/image_analysis/bragg_integration/ProfileIntegrate2D.cpp @@ -39,6 +39,8 @@ Rough BoxSum(const T *image, size_t xpixel, size_t ypixel, int64_t special, int6 int64_t I_sum = 0, n_inner = 0, n_inner_valid = 0; double bkg_sum = 0.0; int n_bkg = 0; + std::vector bkg_vals; + bkg_vals.reserve(128); for (int64_t y = y0; y <= y1; ++y) for (int64_t x = x0; x <= x1; ++x) { const double d2 = (x - r.predicted_x) * (x - r.predicted_x) + (y - r.predicted_y) * (y - r.predicted_y); @@ -52,11 +54,23 @@ Rough BoxSum(const T *image, size_t xpixel, size_t ypixel, int64_t special, int6 if (px == special || px == special + 1 || px == saturation || px == saturation - 1) continue; bkg_sum += static_cast(px); ++n_bkg; + bkg_vals.push_back(static_cast(px)); } } if (n_inner_valid != n_inner || n_bkg <= 5) return out; out.bkg = bkg_sum / n_bkg; + // One high-outlier sigma-clip pass on the background ring. A bandwidth-streaked high-resolution spot + // (or a neighbour) leaks into the annulus and biases the mean high, which over-subtracts and drives + // the merged high-resolution intensities negative; rejecting pixels above mean + 3*sqrt(mean) removes + // that contamination (and barely touches a clean Poisson background, where ~0.1% exceed the cut). + { + const double thr = out.bkg + 3.0 * std::sqrt(std::max(out.bkg, 1.0)); + double s = 0.0; + int n = 0; + for (const double v : bkg_vals) if (v <= thr) { s += v; ++n; } + if (n > 5) out.bkg = s / n; + } out.I = static_cast(I_sum) - static_cast(n_inner) * out.bkg; out.sigma = std::max(1.0, out.I * min_sigma_ratio); if (I_sum > 0)