DetectorModules: panel look-up by id, plane -> panel, plane bounds and a per-pixel walk; mask size checked

Stage 0 of analysing panel detectors in their stored layout. The table of modules already is the panel
table; what the consumers converted in the next stages need from it:

- ImageToPlane / PlaneToImage on a known panel, no search; panel -1 (empty table) is the identity.
- PlaneToPanel: the panel holding a plane point and its image position, nullopt in a gap (never passed
  through, unlike the existing PlaneToImage).
- PlaneBounds: the box of all panel footprints on the plane; with an empty table exactly the image
  rectangle today's prediction bounds test uses.
- ForEachPixel: f(panel, index, px, py) panel by panel, row by row, gaps not visited; with an empty
  table the plain row-major loop with plane == image.

The image preprocessors (CPU, GPU) and the Bragg engines (CPU, GPU; Metal derives from CPU) now refuse
a pixel mask whose size is not the image's, instead of reading a mask of another size out of range.
Every reader sizes the experiment's detector, and so the mask, by the image it reads, so nothing that
runs today is refused. Nothing calls the new functions yet: no output changes.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_015CeEGgpKdfkCTJb3YNMJcu
This commit is contained in:
leonarski_fandClaude Opus 5.5 committed 2026-10-21 09:41:35 +02:00
1 parent b9ada7914b
commit a8202789f5
6 files changed
+104 -2

No files matched your search

+41
View File
@@ -112,6 +112,47 @@ std::pair<double, double> DetectorModules::PlaneToImage(double px, double py) co
return modules[m].PlaneToImage(px, py);
}
std::pair<double, double> DetectorModules::ImageToPlane(int64_t panel, double x, double y) const {
if (panel < 0)
return {x, y};
return modules.at(panel).ImageToPlane(x, y);
}
std::pair<double, double> DetectorModules::PlaneToImage(int64_t panel, double px, double py) const {
if (panel < 0)
return {px, py};
return modules.at(panel).PlaneToImage(px, py);
}
std::optional<PanelPoint> DetectorModules::PlaneToPanel(double px, double py) const {
if (modules.empty())
return PanelPoint{-1, px, py};
for (size_t i = 0; i < modules.size(); i++) {
const auto [x, y] = modules[i].PlaneToImage(px, py);
if (modules[i].ContainsImagePixel(x, y))
return PanelPoint{static_cast<int64_t>(i), x, y};
}
return std::nullopt;
}
PlaneBox DetectorModules::PlaneBounds(int64_t width, int64_t height) const {
if (modules.empty())
return {0.0, 0.0, static_cast<double>(width), static_cast<double>(height)};
PlaneBox box{INFINITY, INFINITY, -INFINITY, -INFINITY};
for (const auto &m: modules) {
const double xa = static_cast<double>(m.x0) - 0.5, xb = static_cast<double>(m.x0 + m.width) - 0.5;
const double ya = static_cast<double>(m.y0) - 0.5, yb = static_cast<double>(m.y0 + m.height) - 0.5;
for (const auto &[x, y]: {std::pair{xa, ya}, std::pair{xb, ya}, std::pair{xa, yb}, std::pair{xb, yb}}) {
const auto [px, py] = m.ImageToPlane(x, y);
box.x0 = std::min(box.x0, px);
box.y0 = std::min(box.y0, py);
box.x1 = std::max(box.x1, px);
box.y1 = std::max(box.y1, py);
}
}
return box;
}
void DetectorModules::Apply(const std::vector<ModuleShift> &shifts) {
for (auto &m: modules) {
const auto [cx, cy] = m.PlaneCentre();
+47
View File
@@ -5,6 +5,7 @@
#include <cstddef>
#include <cstdint>
#include <optional>
#include <utility>
#include <vector>
@@ -39,7 +40,20 @@ struct DetectorModule {
void Apply(const ModuleShift &shift);
};
// A point in the image rugnux reads, on panel `panel` (an index into DetectorModules; -1 when the table is
// empty and the whole image is one flat panel).
struct PanelPoint {
int64_t panel;
double x, y;
};
// An axis-aligned box on the detector plane, x0 <= x < x1, y0 <= y < y1.
struct PlaneBox {
double x0, y0, x1, y1;
};
// The modules of one detector. Empty means the image is one flat detector: every mapping is the identity.
// A module is a panel: one rigid box of the image with one placement. Panel ids are indices into the table.
class DetectorModules {
std::vector<DetectorModule> modules;
public:
@@ -64,6 +78,39 @@ public:
// Plane -> image; a point on no placed module is left as it is.
[[nodiscard]] std::pair<double, double> PlaneToImage(double px, double py) const;
// Image <-> plane on a known panel, no search. Panel -1 (the empty table) is the identity.
[[nodiscard]] std::pair<double, double> ImageToPlane(int64_t panel, double x, double y) const;
[[nodiscard]] std::pair<double, double> PlaneToImage(int64_t panel, double px, double py) const;
// The panel whose placed footprint holds a plane point, and where that point is in the image; nullopt in
// a gap or off every panel. Empty table: panel -1 and the point itself, wherever it is - the caller
// bounds it by the image, as without panels.
[[nodiscard]] std::optional<PanelPoint> PlaneToPanel(double px, double py) const;
// Box on the plane enclosing every panel's pixels (their footprints, half a pixel beyond the edge pixel
// centres): a quick reject before PlaneToPanel. Empty table: [0, width) x [0, height), the image
// rectangle the bounds tests without panels use.
[[nodiscard]] PlaneBox PlaneBounds(int64_t width, int64_t height) const;
// Calls f(panel, index, px, py) for every image pixel on a panel - panel by panel, row by row - with
// index = y * width + x in the image and (px, py) its centre on the plane. Pixels on no panel (gaps) are
// not visited. Empty table: one panel -1 that is the whole image, plane == image.
template<class F>
void ForEachPixel(int64_t width, int64_t height, F &&f) const {
if (modules.empty()) {
for (int64_t y = 0; y < height; y++)
for (int64_t x = 0; x < width; x++)
f(int64_t(-1), static_cast<size_t>(y * width + x), static_cast<double>(x), static_cast<double>(y));
return;
}
for (size_t p = 0; p < modules.size(); p++) {
const DetectorModule &m = modules[p];
for (int64_t y = m.y0; y < m.y0 + m.height; y++)
for (int64_t x = m.x0; x < m.x0 + m.width; x++) {
const auto [px, py] = m.ImageToPlane(static_cast<double>(x), static_cast<double>(y));
f(static_cast<int64_t>(p), static_cast<size_t>(y * width + x), px, py);
}
}
}
// Each module moved by the shift whose centre is nearest its own plane centre (and within it).
void Apply(const std::vector<ModuleShift> &shifts);
// The modules as refinement starting points: their plane centres, no move.
@@ -78,7 +78,11 @@ BraggIntegrationEngineCPU::BraggIntegrationEngineCPU(const DiffractionExperiment
: BraggIntegrationEngine(experiment),
refl_mask(static_cast<int>(xpixel), static_cast<int>(ypixel), 0),
owner(static_cast<int>(xpixel), static_cast<int>(ypixel), BRAGG_OWNER_NONE),
static_mask(mask.GetPackedMask()) {}
static_mask(mask.GetPackedMask()) {
if (mask.GetMask().size() != npixel)
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
"Pixel mask size doesn't match the image");
}
void BraggIntegrationEngineCPU::MarkNeighbours(const std::vector<Reflection> &predicted,
const std::vector<uint32_t> &order) {
@@ -834,6 +834,9 @@ BraggIntegrationEngineGPU::BraggIntegrationEngineGPU(const DiffractionExperiment
d_global_n(1),
d_invd2(2),
d_counts(COUNT_SLOTS) {
if (mask.GetMask().size() != npixel)
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
"Pixel mask size doesn't match the image");
threads = 128;
// On the engine's stream (the member: the parameter was moved from), where the kernels counting
// into it run - the NULL stream is not ordered before those.
@@ -9,7 +9,11 @@
ImagePreprocessorCPU::ImagePreprocessorCPU(const DiffractionExperiment &experiment, const PixelMask &mask)
: ImagePreprocessor(experiment),
mask_bits(mask.GetPackedMask()) {}
mask_bits(mask.GetPackedMask()) {
if (mask.GetMask().size() != npixels)
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
"Pixel mask size doesn't match the image");
}
ImageStatistics ImagePreprocessorCPU::Analyze(ImagePreprocessorBuffer &processed_image, const uint8_t *image_ptr, CompressedImageMode image_mode) {
if (processed_image.size() != npixels)
@@ -329,6 +329,9 @@ ImagePreprocessorGPU::ImagePreprocessorGPU(const DiffractionExperiment &experime
copy_image_to_host(copy_image_to_host),
gpu_stats(1),
cpu_stats(1) {
if (mask.GetMask().size() != npixels)
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
"Pixel mask size doesn't match the image");
// Setup mask. The byte-per-pixel form and its checksum come from the PixelMask, which derives them
// whenever the mask changes: they are the same for every worker, and deriving them walks every
// pixel of the detector - 18 million of them on a 16 Mpx one, per engine, with an engine built per