Files
leonarski_fandClaude Opus 5.5 a609f6ba1c rugnux two-pass: the integration radius decided on the pre-pass frames, before the canonical pass
Where the measured spot width widened the integration radius, the canonical pass
integrated the whole sweep at the widened radius, and only then did the starvation
guard read how many background rings the neighbouring predictions crowd; over the
bound (1.13%) the pass was thrown away and run again at the fixed radius.

The geometry pre-pass now counts that same quantity at the widened radii without
integrating there (BraggIntegrationEngineCPU::CountNeighbourCrowding: the neighbour
mask and the ring pixels the run's mask leaves readable, no image) on 20 frames
spread over the pre-pass, with the standard error of the ratio over those frames.
Over the bound, the canonical pass integrates at the fixed radius from the outset.
A count within NEAR_TIE_SIGMA of the bound is a near tie in the part-B sense: the
canonical pass counts the widened radii again on 20 frames over the whole sweep,
and where that puts them under the bound the run is made again deciding late. A
run deciding late counts nothing early and decides as before. Under the bound the
canonical pass integrates at the widened radius and the measured guard stays as it
was, so a pre-pass count that is too low (its prediction holds fewer tails than
the canonical pass's, or it stood on another lattice class) costs only the pass
it cost before. No new threshold.

Battery against 46be3f647 (the 14 open sets that re-ran at the fixed radius, the
9 open/in-house sets closest under the bound, the ci tier, 2 private sets):
- p.hkl and MTZ data identical on all 51 open/in-house sets and on one private set.
- The 14 re-running sets: rugnux wall 534 -> 447 s (7qij -27, 4nwv -12, 7yzx -10,
  9fhc -9, 5src -6 s, on a shared machine). 8sqq and 8xtf count under the bound on
  the pre-pass and re-run as before. 4nwv (1.69 +- 0.28%) and one other set were
  near ties, re-checked on the whole sweep (1.72%, 1.61%), ending where they did.
- One private, densely crowded set: same space group and d_min, R_meas 28.6 ->
  39.3%. Its old re-run re-indexed starting from the spot budget the thrown-away
  pass had measured (152 of 1000). With that budget given (--max-spots 152), the
  new run gives 30.6%. The difference comes from that path, not from the radius.
- The count costs ~13 CPU-s on the heaviest set (65k predictions a frame). Over
  the 22 ci sets that count, pre-pass integration took 45.2 -> 45.8 s.
- A build forcing the near tie and its overturn on 4nwv re-checks on the whole
  sweep (1.80% against the engine's 1.81%) and makes the run again deciding late.
  It ends with the same p.hkl.

Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01EizBKhTqrAYago9KJG3dA3
2026-10-11 13:01:23 +02:00

130 lines
7.0 KiB
C++

// SPDX-FileCopyrightText: 2024 Filip Leonarski, Paul Scherrer Institute <filip.leonarski@psi.ch>
// SPDX-License-Identifier: GPL-3.0-only
#pragma once
#include <mutex>
#include "../common/JFJochMessages.h"
#include "../common/DiffractionExperiment.h"
#include "../common/AzimuthalIntegrationMapping.h"
#include "../common/PixelMask.h"
#include "../common/AzimuthalIntegrationProfile.h"
#include "bragg_prediction/BraggPrediction.h"
#include "bragg_integration/BraggIntegrationEngine.h"
#include "bragg_integration/BraggIntegrationEngineCPU.h"
#include "spot_finding/ImageSpotFinder.h"
#include "spot_finding/AdaptiveSpotFinderCPU.h"
#include "indexing/IndexerThreadPool.h"
#include "azint/AzIntEngine.h"
#include "roi/ROIIntegration.h"
#include "IndexAndRefine.h"
#include "image_preprocessing/ImagePreprocessor.h"
#include "image_preprocessing/ImagePreprocessorBuffer.h"
#include "image_preprocessing/ImagePreprocessorCPU.h"
class CudaStream;
class AdaptiveSpotFinderGPU;
// MXAnalysisWithoutFPGA is not thread safe - it has to owned by a single thread
class MXAnalysisWithoutFPGA {
const DiffractionExperiment &experiment;
const AzimuthalIntegrationMapping &integration;
std::vector<uint8_t> decompression_buffer;
std::unique_ptr<ImagePreprocessor> preprocessor;
// The preprocessor, where it is the CPU one.
ImagePreprocessorCPU *preprocessor_cpu = nullptr;
size_t npixels;
size_t xpixels;
// Built on first use: the fused adaptive finder produces the azimuthal profile as a by-product,
// so on the rugnux path this engine is constructed and then never run.
std::unique_ptr<AzIntEngine> azint;
AzIntEngine &AzInt();
std::unique_ptr<ROIIntegration> roi;
// Built on first use. Which finder an image takes arrives with its SpotFindingSettings, and
// with adaptive detection on - the default everywhere but the broker - this one is never asked
// for; on the GPU it is ~14 MB and 15 device allocations per worker.
std::unique_ptr<ImageSpotFinder> spotFinder;
ImageSpotFinder &FixedThresholdFinder();
// Self-calibrating finder, used when spot settings request adaptive detection. Kept alongside the
// default finder because the choice arrives with the per-image settings, not at construction. It is
// an AdaptiveSpotFinderCPU by default; on the GPU path, when the fused engine is enabled (rugnux
// offline only), it is instead an AdaptiveSpotFinderGPU that also computes the azimuthal profile,
// aliased through fused_adaptive so Analyze() can take that profile and skip the separate azint pass.
std::unique_ptr<ImageSpotFinder> adaptiveSpotFinder;
AdaptiveSpotFinderGPU *fused_adaptive = nullptr;
// The CPU finder, where it gives the profile too (the azimuthal integration being on the CPU).
AdaptiveSpotFinderCPU *fused_adaptive_cpu = nullptr;
const bool enable_fused_adaptive_gpu;
IndexAndRefine &indexer;
std::unique_ptr<BraggPrediction> prediction;
std::unique_ptr<BraggIntegrationEngine> bragg_engine;
// What the supercell probe's integrations added to bragg_engine's counts (see BraggCounts).
BraggIntegrationCounts probe_counts;
// Counts the neighbour crowding of the integrated reflections at another stencil (CountCrowdingAt),
// on the images for which count_crowding is on.
std::unique_ptr<BraggIntegrationEngineCPU> crowding_engine;
bool count_crowding = false;
std::unique_ptr<ImagePreprocessorBuffer> preprocessor_buffer;
const PixelMask &mask;
// Decompress the image into decompression_buffer (or read it straight from the message, when it is
// not compressed) and return where it landed.
const uint8_t *Decompress(const CompressedImage &image);
// The CPU preprocessing. A bitshuffled image is decoded a block (~32 KB) at a time and each block
// is preprocessed while it is in cache, then - with ring_pass - put through the adaptive finder's
// plain ring pass as well, so the image is never written out whole before it is preprocessed and
// the preprocessed pixels are read back from cache, not from memory. The blocks come in pixel
// order, so every sum is taken in the same order as by separate passes.
ImageStatistics PreprocessCPU(const CompressedImage &image, bool ring_pass);
// Pixels outside the resolution limits, bit-packed. Built by the integration mapping, which is
// shared by every worker's engine and hands out the same mask to all of them.
std::shared_ptr<const std::vector<uint32_t>> mask_resolution;
// The limits mask_resolution was built for. Kept as the OPTIONAL the caller passed, so an unset
// high-resolution limit compares equal to itself and the mask is not rebuilt on every image.
std::optional<float> mask_high_res;
std::optional<float> mask_low_res;
void UpdateMaskResolution(const SpotFindingSettings& settings);
#ifdef JFJOCH_USE_CUDA
std::shared_ptr<CudaStream> stream; // kept so RebuildROI() can recreate the GPU ROI engine
#endif
public:
// enable_fused_adaptive_gpu turns on the fused GPU azint+adaptive spot finder (only takes effect on
// the GPU path with adaptive detection). The rugnux offline path and the interactive viewer enable
// it by default, as does the online receiver. It only changes performance - the fused engine
// reproduces the CPU finder's spots. Note it also decides whether the preprocessed image is copied
// back to the host each frame: that copy exists only for a CPU engine to read, and with the flag on
// no CPU engine is built, so the copy is skipped.
MXAnalysisWithoutFPGA(const DiffractionExperiment &experiment, const AzimuthalIntegrationMapping &integration,
const PixelMask &mask, IndexAndRefine &indexer, bool enable_fused_adaptive_gpu = false);
void Analyze(DataMessage &output, AzimuthalIntegrationProfile &profile, const SpotFindingSettings &spot_finding_settings);
// Surgical ROI-only paths used when a full re-analysis is not wanted: rebuild the
// ROI engine after the ROI set changes, recompute ROIs after preprocessing a new
// image (reanalyze off), or just rerun ROIs on the current preprocessed image (an
// interactive ROI move). A full Analyze() already computes ROIs, so needs nothing.
void RebuildROI();
void AnalyzeROIOnly(DataMessage &output);
void RunROIOnly(DataMessage &output);
// What this worker's Bragg integrator counted (BraggIntegrationCounts). Each worker builds its own
// analysis, so a caller that wants the run's totals sums this over the workers it started.
[[nodiscard]] BraggIntegrationCounts BraggCounts() const;
// Also count, on the images Analyze is given while CountCrowding(true), how crowded the background
// rings of the integrated reflections would be at `settings` rather than the experiment's own
// (BraggIntegrationEngineCPU::CountNeighbourCrowding), read back with CrowdingCounts.
void CountCrowdingAt(const BraggIntegrationSettings &settings);
void CountCrowding(bool on) { count_crowding = on; }
[[nodiscard]] BraggIntegrationCounts CrowdingCounts() const;
};