diff --git a/image_analysis/hot_pixels/HotPixels.cpp b/image_analysis/hot_pixels/HotPixels.cpp index 29937caec..79736c80b 100644 --- a/image_analysis/hot_pixels/HotPixels.cpp +++ b/image_analysis/hot_pixels/HotPixels.cpp @@ -203,7 +203,8 @@ void HotPixelFinder::AddImage(const int32_t *image, std::vector &scratc for (size_t i = begin; i < end; i++) { const int32_t k = key[i], v = image[i]; if (k < 0) continue; - if (v == INT32_MIN) { + // The error value and the saturation limit alike: no count (see the class comment). + if (v == INT32_MIN || v == INT32_MAX) { n_error[i]++; if (ring_ok[k / SECTORS]) { // counted in its ring-sector's frames, but not valid here n_error_ring_ok[i]++; @@ -213,7 +214,7 @@ void HotPixelFinder::AddImage(const int32_t *image, std::vector &scratc } if (!ring_ok[k / SECTORS]) continue; sum_value[i] += v; - if (v == INT32_MAX || static_cast(v) > threshold[k]) + if (static_cast(v) > threshold[k]) n_lit[i]++; } } diff --git a/image_analysis/hot_pixels/HotPixels.h b/image_analysis/hot_pixels/HotPixels.h index e364871e0..7320979bc 100644 --- a/image_analysis/hot_pixels/HotPixels.h +++ b/image_analysis/hot_pixels/HotPixels.h @@ -31,8 +31,10 @@ // is above the Poisson bound: a weaker one cannot make an outlier. Alone or as part of a patch - a // patch lit through the whole sweep is stationary in the lab, so it is not a reflection of the // rotating crystal, and a reflection crossing it reads the patch. A pixel holding the -// detector's error value on most frames is masked with it: every frame already treats it as invalid, -// but the static mask did not know. +// detector's error value or its saturation limit on most frames is masked with it: every frame already +// treats it as invalid, but the static mask did not know. A pixel stuck at the limit is no overload - +// a CCD can have thousands, at 65535 on every frame - and unmasked, every reflection whose disk it +// falls in is flagged as overloaded. // // A DEAD pixel is the mirror of a hot one: it reads on average at most half of its own sector's level, // and its mean deficit is above the bound a hot pixel's excess has to clear. Nothing but a defect keeps @@ -85,7 +87,7 @@ public: struct Result { std::vector mask; // non-zero = masked (1 hot, 2 error value, 3 dead), converted geometry size_t hot = 0; // persistent above the ring - size_t error = 0; // holding the error value on most frames + size_t error = 0; // holding the error value or the saturation limit on most frames size_t dead = 0; // reading at most half of its sector's level uint32_t frames = 0; }; diff --git a/image_analysis/hot_pixels/HotPixelsGPU.cu b/image_analysis/hot_pixels/HotPixelsGPU.cu index a5725533b..61dbc3474 100644 --- a/image_analysis/hot_pixels/HotPixelsGPU.cu +++ b/image_analysis/hot_pixels/HotPixelsGPU.cu @@ -132,7 +132,7 @@ __global__ void accumulate_kernel(const int32_t *__restrict__ image, const int32 i += static_cast(blockDim.x) * gridDim.x) { const int32_t k = key[i], v = image[i]; if (k < 0) continue; - if (v == INT32_MIN) { + if (v == INT32_MIN || v == INT32_MAX) { // the error value or the saturation limit: no count n_error[i]++; if (ring_ok[k / sectors]) { n_error_ring_ok[i]++; @@ -142,7 +142,7 @@ __global__ void accumulate_kernel(const int32_t *__restrict__ image, const int32 } if (!ring_ok[k / sectors]) continue; sum_value[i] += v; - if (v == INT32_MAX || static_cast(v) > threshold[k]) + if (static_cast(v) > threshold[k]) n_lit[i]++; } } diff --git a/rugnux/RugnuxPreScan.cpp b/rugnux/RugnuxPreScan.cpp index 9603581df..6786be5be 100644 --- a/rugnux/RugnuxPreScan.cpp +++ b/rugnux/RugnuxPreScan.cpp @@ -1218,7 +1218,7 @@ void Rugnux::MaskDefectivePixels(int start_image, const std::vector &sample pixel_mask_.LoadHotPixelMask(experiment_, mask); const size_t masked = defects->hot + defects->error + defects->dead; logger.Info("Defective pixels: {} masked from {} pre-scan frames - {} lit above their ring on more " - "frames than one reflection or chance explains, {} holding the detector's error value, {} dead", + "frames than one reflection or chance explains, {} holding the detector's error value or saturation limit, {} dead", masked, defects->frames, defects->hot, defects->error, defects->dead); Note(fmt::format("{} defective pixel{} masked", masked, masked == 1 ? "" : "s")); }; diff --git a/tests/HotPixelFinderTest.cpp b/tests/HotPixelFinderTest.cpp index a4015e8b8..da00f00cd 100644 --- a/tests/HotPixelFinderTest.cpp +++ b/tests/HotPixelFinderTest.cpp @@ -25,9 +25,10 @@ namespace { // A Poisson background with a powder ring, and on it four kinds of pixel: one that reads 60 counts // high on every frame, one only 20 high - persistent, but too weak to make an outlier - one that holds -// the error value on every frame, and one crossed by a genuine reflection. The last sits close to the rotation axis, where one reflection stays on a pixel for a +// the error value on every frame, one stuck at the saturation limit, and one crossed by a genuine +// reflection, saturated on a few frames. The last sits close to the rotation axis, where one reflection stays on a pixel for a // long stretch of rotation - here nine consecutive sampled frames, more than the chance bound allows -// but fewer than the one-reflection bound at that zeta. Only the first and the third are masked - and +// but fewer than the one-reflection bound at that zeta. Only the first, third and fourth are masked - and // a 3x3 patch reading high on every frame, which is stationary in the lab and so no reflection either. // Three dark pixels besides: one dead on the powder ring, which is masked, one on the ring reading 70% // of it, and one dead on the 5-count background, whose deficit is too small to make an outlier. @@ -40,7 +41,7 @@ TEST_CASE("HotPixelFinder_PersistentPixelNotBragg", "[HotPixelFinder]") { HotPixelFinder finder(x, pixel_mask, 4); constexpr int HOT_X = 200, HOT_Y = 200, WARM_X = 190, WARM_Y = 60, ERR_X = 60, ERR_Y = 190; - constexpr int BRAGG_X = 40, BRAGG_Y = 115; + constexpr int BRAGG_X = 40, BRAGG_Y = 115, STUCK_X = 100, STUCK_Y = 30; constexpr int PATCH_X = 150, PATCH_Y = 225; constexpr int DEAD_X = C + 82, DEAD_Y = C, DIM_X = C, DIM_Y = C + 82, DEAD_BKG_X = 20, DEAD_BKG_Y = 20; std::mt19937 rng(1); @@ -57,11 +58,12 @@ TEST_CASE("HotPixelFinder_PersistentPixelNotBragg", "[HotPixelFinder]") { frame[I(HOT_X, HOT_Y)] += 60; frame[I(WARM_X, WARM_Y)] += 20; frame[I(ERR_X, ERR_Y)] = INT32_MIN; + frame[I(STUCK_X, STUCK_Y)] = INT32_MAX; for (int dy = -1; dy <= 1; dy++) for (int dx = -1; dx <= 1; dx++) frame[I(PATCH_X + dx, PATCH_Y + dy)] += 60; if (f >= 20 && f < 29) - frame[I(BRAGG_X, BRAGG_Y)] += 500; + frame[I(BRAGG_X, BRAGG_Y)] = f < 23 ? INT32_MAX : frame[I(BRAGG_X, BRAGG_Y)] + 500; frame[I(DEAD_X, DEAD_Y)] = 0; frame[I(DIM_X, DIM_Y)] = static_cast(0.7 * frame[I(DIM_X, DIM_Y)]); frame[I(DEAD_BKG_X, DEAD_BKG_Y)] = 0; @@ -72,6 +74,7 @@ TEST_CASE("HotPixelFinder_PersistentPixelNotBragg", "[HotPixelFinder]") { CHECK(result.frames == NFRAMES); CHECK(result.mask[I(HOT_X, HOT_Y)] == 1); // 1 = hot, 2 = error value, 3 = dead CHECK(result.mask[I(ERR_X, ERR_Y)] == 2); + CHECK(result.mask[I(STUCK_X, STUCK_Y)] == 2); CHECK(result.mask[I(WARM_X, WARM_Y)] == 0); CHECK(result.mask[I(BRAGG_X, BRAGG_Y)] == 0); for (int dy = -1; dy <= 1; dy++) @@ -81,7 +84,7 @@ TEST_CASE("HotPixelFinder_PersistentPixelNotBragg", "[HotPixelFinder]") { CHECK(result.mask[I(DIM_X, DIM_Y)] == 0); CHECK(result.mask[I(DEAD_BKG_X, DEAD_BKG_Y)] == 0); CHECK(result.hot == 10); - CHECK(result.error == 1); + CHECK(result.error == 2); CHECK(result.dead == 1); }