From 2041ce2bd5ede44b345854228ca081f5d9c9c8b9 Mon Sep 17 00:00:00 2001
From: leonarski_f <filip.leonarski@psi.ch>
Date: Wed, 24 Jun 2026 22:42:25 +0200
Subject: [PATCH] Reject the +/-1 masked/saturated sentinel band in 2D
 integration

Masked/error pixels carry the int type minimum and saturated the maximum,
but the lossy codec can nudge them inward by one (masked observed as
INT32_MIN+1 in the decompressed data). The integrators only checked the
exact extremes, so a shifted sentinel in a reflection's background ring was
treated as a valid pixel -> garbage background and intensity for any
reflection whose box clips a module gap. Reject the +/-1 band too (real
calibrated counts never approach the type extremes). Neutral on the
well-centred lyso test set; a correctness fix for gap-clipping reflections.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
---
 image_analysis/bragg_integration/BraggIntegrate2D.cpp   | 5 ++++-
 image_analysis/bragg_integration/ProfileIntegrate2D.cpp | 8 ++++----
 2 files changed, 8 insertions(+), 5 deletions(-)

diff --git a/image_analysis/bragg_integration/BraggIntegrate2D.cpp b/image_analysis/bragg_integration/BraggIntegrate2D.cpp
index 22339929..d8ea5dce 100644
--- a/image_analysis/bragg_integration/BraggIntegrate2D.cpp
+++ b/image_analysis/bragg_integration/BraggIntegrate2D.cpp
@@ -88,7 +88,10 @@ void IntegrateReflection(Reflection &r, const T *image, const std::vector<uint8_
             if (dist_sq < r_1_sq)
                 I_npixel_inner++;
 
-            if (pixel == special_value || pixel == saturation)
+            // Masked/error pixels carry the type minimum and saturated the type maximum, but the lossy
+            // codec can nudge them inward by one (masked observed as INT32_MIN+1), so reject the +/-1
+            // band too - real calibrated counts never get anywhere near the type extremes.
+            if (pixel == special_value || pixel == special_value + 1 || pixel == saturation || pixel == saturation - 1)
                 continue;
 
             if (dist_sq < r_1_sq) {
diff --git a/image_analysis/bragg_integration/ProfileIntegrate2D.cpp b/image_analysis/bragg_integration/ProfileIntegrate2D.cpp
index 9d314d22..c43a5d2d 100644
--- a/image_analysis/bragg_integration/ProfileIntegrate2D.cpp
+++ b/image_analysis/bragg_integration/ProfileIntegrate2D.cpp
@@ -45,11 +45,11 @@ Rough BoxSum(const T *image, size_t xpixel, size_t ypixel, int64_t special, int6
             const auto px = image[y * xpixel + x];
             if (d2 < r1_sq) {
                 ++n_inner;
-                if (px == special || px == saturation) continue;
+                if (px == special || px == special + 1 || px == saturation || px == saturation - 1) continue;
                 I_sum += px;
                 ++n_inner_valid;
             } else if (d2 >= r2_sq && d2 < r3_sq) {
-                if (px == special || px == saturation) continue;
+                if (px == special || px == special + 1 || px == saturation || px == saturation - 1) continue;
                 bkg_sum += static_cast<double>(px);
                 ++n_bkg;
             }
@@ -124,7 +124,7 @@ std::vector<Reflection> ProfileIntegrateInternal(const DiffractionExperiment &ex
                 const int64_t x = rh.cx + dx, y = rh.cy + dy;
                 if (x < 0 || y < 0 || x >= static_cast<int64_t>(xpixel) || y >= static_cast<int64_t>(ypixel)) continue;
                 const auto px = ptr[y * xpixel + x];
-                if (px == special || px == saturation) continue;
+                if (px == special || px == special + 1 || px == saturation || px == saturation - 1) continue;
                 const double v = (static_cast<double>(px) - rh.bkg) / rh.I;
                 shell_grid[rh.shell][idx(dx, dy)] += v;
                 global_grid[idx(dx, dy)] += v;
@@ -191,7 +191,7 @@ std::vector<Reflection> ProfileIntegrateInternal(const DiffractionExperiment &ex
                     const int64_t x = rh.cx + dx, y = rh.cy + dy;
                     if (x < 0 || y < 0 || x >= static_cast<int64_t>(xpixel) || y >= static_cast<int64_t>(ypixel)) continue;
                     const auto px = ptr[y * xpixel + x];
-                    if (px == special || px == saturation) continue;
+                    if (px == special || px == special + 1 || px == saturation || px == saturation - 1) continue;
                     const double v = B + std::max(0.0, I) * Pp;
                     num += Pp * (static_cast<double>(px) - rh.bkg) / v;
                     den += Pp * Pp / v;