Live stills: the geometry is fixed when the pre-pass ends

The running geometry refinement moved the distance during the measurement, so the radial profile's
binning and the frames merged depended on when the fit landed. It is removed, with the analysis threads'
rebuild and the frames held back from the merge: every frame is analysed and merged at the pre-pass's
geometry. Where the pre-pass could not refine the geometry, the header's distance stays for the whole
run, with a warning in the log and a note in the live status. The offline pipeline is not touched.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
This commit is contained in:
leonarski_fandClaude Opus 5.5 committed 2026-11-06 16:48:23 +01:00
1 parent df861de7fb
commit d07e6e2aa9
9 files changed
+37 -351

No files matched your search

+5 -9
View File
@@ -266,15 +266,11 @@ all: on fewer, too few index for the geometry refinement to run, and the stated
buffer has to hold them, and the frames that arrive while the pre-pass runs, in half its slots (for a JUNGFRAU
8M, a slot is the size of an uncompressed frame); a smaller buffer starts the pre-pass on as many as it holds.
Even 500 frames give the pre-pass only a few dozen indexed ones at a serial hit rate, and none where the
first frames are empty (a fixed-target chip whose crystals sit further in). So, where a cell is given and the
pre-pass refined the geometry on fewer than 200 frames, the centre and distance are refined again on the frames
as they are analysed, from their spots and lattices, without holding their images: once 200 strong frames (20
indexed spots or more) are in, or 5000 indexed, the geometry is fitted, as rugnux fits it on its sample, and the
analysis switches to it; a second round does the same on frames analysed at that geometry. Two fits in all, a
fixed cost. The frames analysed before the last switch are written as analysed but left out of the merge (their
reflections were predicted at the earlier geometry, and the images are no longer there to integrate again);
the switches are notes of the live status, and the end message carries the geometry the run ended at.
The geometry is fixed when the pre-pass ends and does not change during the collection, so every frame is analysed
and merged at it and the radial profile has one binning throughout. Where the pre-pass could not refine it (too few
of its frames indexed, as with a small image buffer, a low hit rate or a chip whose first frames are empty, or no unit
cell given), the distance in the header is used for the whole run: this is a warning in the log and a note of the
live status, and a wrong header distance then costs resolution and indexed frames.
Frames arriving during the pre-pass, or faster than the analysis, wait in the image buffer, in at most half of
its slots; a frame that finds no free slot is dropped and counted as for a full image buffer, so analysis
-6
View File
@@ -1128,12 +1128,6 @@ IntegrationOutcome IndexAndRefine::TakeIntegrationOutcome(int64_t image) {
return std::move(integration_outcome.at(image));
}
void IndexAndRefine::DropIntegrationOutcome(int64_t image) {
const std::unique_lock ul(reflections_mutex);
integration_outcome.at(image) = {};
unit_cells.at(image).reset();
}
std::vector<IntegrationOutcome> &IndexAndRefine::GetIntegrationOutcome() {
return integration_outcome;
}
-2
View File
@@ -309,7 +309,5 @@ public:
// The retained outcome of one image, moved out: what is left in its place is empty. Safe beside the
// analysis of other images.
IntegrationOutcome TakeIntegrationOutcome(int64_t image);
// The retained outcome of one image dropped, and its cell with it: it is not in the consensus.
void DropIntegrationOutcome(int64_t image);
const std::vector<IntegrationOutcome> &GetIntegrationOutcome() const;
};
+28 -126
View File
@@ -207,15 +207,15 @@ void StillsGPUAnalysis::PrePassThread() {
setup.mapping = std::make_unique<AzimuthalIntegrationMapping>(setup.experiment, setup.mask);
// Retained for the merge, which takes each image's reflections as soon as it is analysed.
indexer = std::make_unique<IndexAndRefine>(setup.experiment, indexing_thread_pool, merge_frames, true);
// The pre-pass's fit rests on the frames of its sample that indexed; where those are fewer than the
// running refinement bundles, the geometry is refined again, from where the pre-pass left it, on the
// frames as they are analysed.
if (geometry_frames < StillsRunningGeometry::STRONG_FRAMES && setup.experiment.GetUnitCell()) {
running_geometry = std::make_unique<StillsRunningGeometry>(setup.experiment, MERGE_THREADS);
frame_epoch.assign(experiment.GetImageNum(), 0);
geometry_state = GEOMETRY_PENDING;
logger.Info("The pre-pass refined the geometry on {} frames: it is refined again on {} strong frames as "
"they are analysed", geometry_frames, StillsRunningGeometry::STRONG_FRAMES);
if (geometry_frames == 0) {
const std::string text = fmt::format(
"The pre-pass could not refine the geometry (too few frames indexed, or no unit cell given): the "
"whole run is analysed at the distance in the header, {:.2f} mm, beam centre ({:.2f}, {:.2f}) px",
setup.experiment.GetDetectorDistance_mm(), setup.experiment.GetBeamX_pxl(),
setup.experiment.GetBeamY_pxl());
logger.Warning("{}", text);
if (status)
status->Note(text);
}
if (merge_frames) {
merge_experiment = setup.experiment;
@@ -267,54 +267,32 @@ void StillsGPUAnalysis::AnalysisThread() {
return;
}
// A thread whose analysis cannot be built still passes its frames on, unanalysed. It is built again at
// the geometry the running refinement switches to.
// A thread whose analysis cannot be built still passes its frames on, unanalysed.
std::unique_ptr<MXAnalysisWithoutFPGA> analysis;
int epoch = -1;
const auto build = [&] {
epoch = geometry_epoch;
analysis.reset();
try {
analysis = std::make_unique<MXAnalysisWithoutFPGA>(setup.experiment, *setup.mapping, setup.mask,
*indexer, /*enable_fused_adaptive_gpu=*/true);
} catch (const std::exception &e) {
AnalysisError(std::string("analysis thread could not start, its frames are written unanalysed - ")
+ e.what());
}
};
try {
analysis = std::make_unique<MXAnalysisWithoutFPGA>(setup.experiment, *setup.mapping, setup.mask,
*indexer, /*enable_fused_adaptive_gpu=*/true);
} catch (const std::exception &e) {
AnalysisError(std::string("analysis thread could not start, its frames are written unanalysed - ")
+ e.what());
}
while (auto held = queue.GetBlocking()) {
DataMessage &message = held->message;
std::optional<AzimuthalIntegrationProfile> profile;
bool switch_geometry = false;
const auto start_time = std::chrono::steady_clock::now();
{
std::unique_lock gate(gate_mutex);
gate_cv.wait(gate, [this] { return !switch_waiting; });
}
{
std::shared_lock sl(geometry_mutex);
if (epoch != geometry_epoch)
build();
profile.emplace(*setup.mapping);
try {
if (analysis)
analysis->Analyze(message, *profile, settings());
} catch (const JFJochException &e) {
// The frame is still written, without its analysis.
logger.ErrorException(e);
} catch (const std::exception &e) {
AnalysisError(fmt::format("image {} written unanalysed - {}", message.number, e.what()));
}
if (geometry_state == GEOMETRY_PENDING) {
frame_epoch.at(message.number) = static_cast<uint8_t>(epoch + 1);
switch_geometry = running_geometry->Add(message);
}
profile.emplace(*setup.mapping);
try {
if (analysis)
analysis->Analyze(message, *profile, settings());
} catch (const JFJochException &e) {
// The frame is still written, without its analysis.
logger.ErrorException(e);
} catch (const std::exception &e) {
AnalysisError(fmt::format("image {} written unanalysed - {}", message.number, e.what()));
}
message.processing_time_s = std::chrono::duration<float>(std::chrono::steady_clock::now()
- start_time).count();
if (switch_geometry)
SwitchGeometry();
try {
sink(message, *profile);
if (merge)
@@ -337,63 +315,6 @@ void StillsGPUAnalysis::AnalysisThread() {
}
}
void MoveStillsGeometry(DiffractionExperiment &experiment, float beam_x_pxl, float beam_y_pxl, float distance_mm) {
AzimuthalIntegrationSettings azint = experiment.GetAzimuthalIntegrationSettings();
azint.QRange_recipA(azint.GetLowQ_recipA(), azint.GetHighQ_recipA());
experiment.ImportAzimuthalIntegrationSettings(azint);
experiment.BeamX_pxl(beam_x_pxl).BeamY_pxl(beam_y_pxl).DetectorDistance_mm(distance_mm);
}
void StillsGPUAnalysis::SwitchGeometry() {
// One fit, in the analysis thread that offered the last frame it needed; the others go on meanwhile.
const auto t0 = std::chrono::steady_clock::now();
const auto r = running_geometry->Refine();
const bool last_round = ++geometry_round == GEOMETRY_ROUNDS;
if (!r) {
logger.Warning("Running geometry refinement, round {}: no fit on {} frames ({} strong); the frames are "
"analysed at the geometry they were", geometry_round, running_geometry->Offered(),
running_geometry->Strong());
geometry_state = geometry_epoch > 0 ? GEOMETRY_SWITCHED : GEOMETRY_KEPT;
return;
}
const std::string text = fmt::format(
"Geometry refined on {} analysed frames ({} strong), round {}: beam centre ({:.2f}, {:.2f}) -> "
"({:.2f}, {:.2f}) px, distance {:.2f} -> {:.2f} mm ({:.1f} s)", running_geometry->Offered(),
running_geometry->Strong(), geometry_round, setup.experiment.GetBeamX_pxl(),
setup.experiment.GetBeamY_pxl(), r->beam_x_px, r->beam_y_px, setup.experiment.GetDetectorDistance_mm(),
r->distance_mm, std::chrono::duration<double>(std::chrono::steady_clock::now() - t0).count());
{
std::unique_lock gate(gate_mutex);
switch_waiting = true;
}
{
std::unique_lock ul(geometry_mutex);
MoveStillsGeometry(setup.experiment, static_cast<float>(r->beam_x_px), static_cast<float>(r->beam_y_px),
static_cast<float>(r->distance_mm));
retired_mappings.push_back(std::move(setup.mapping));
setup.mapping = std::make_unique<AzimuthalIntegrationMapping>(setup.experiment, setup.mask);
++geometry_epoch;
if (last_round)
geometry_state = GEOMETRY_SWITCHED;
else
running_geometry = std::make_unique<StillsRunningGeometry>(setup.experiment, MERGE_THREADS);
}
{
std::unique_lock gate(gate_mutex);
switch_waiting = false;
}
gate_cv.notify_all();
{
std::unique_lock ul(live_mutex);
live.beam_x_pxl = r->beam_x_px;
live.beam_y_pxl = r->beam_y_px;
live.distance_mm = r->distance_mm;
}
logger.Info("{}", text);
if (status)
status->Note(text);
}
std::optional<UnitCell> StillsGPUAnalysis::GetConsensusUnitCell() const {
return indexer->GetConsensusUnitCell();
}
@@ -446,27 +367,8 @@ void StillsGPUAnalysis::MergeThread() {
} else if (next > 0 && merge->Added() % PREPARE_EVERY == 0)
merge->Prepare(MERGE_THREADS);
};
// The frames analysed at a geometry the running refinement may still replace wait until it is decided:
// left out where analysed before its last switch, merged as analysed otherwise.
std::vector<int64_t> held_back;
const auto release = [&] {
for (const int64_t h: held_back)
if (frame_epoch[h] - 1 < geometry_epoch) {
indexer->DropIntegrationOutcome(h);
add(h, IntegrationOutcome{});
} else
add(h, indexer->TakeIntegrationOutcome(h));
held_back.clear();
};
for (int64_t n = merge_queue.GetBlocking(); n >= 0; n = merge_queue.GetBlocking()) {
if (!frame_epoch.empty() && frame_epoch[n] > 0)
held_back.push_back(n);
else
add(n, indexer->TakeIntegrationOutcome(n));
if (geometry_state != GEOMETRY_PENDING && !held_back.empty())
release();
}
release();
for (int64_t n = merge_queue.GetBlocking(); n >= 0; n = merge_queue.GetBlocking())
add(n, indexer->TakeIntegrationOutcome(n));
if (have_cell()) {
merged = merge->Final(static_cast<int>(AvailableCpus()));
logger.Info("Stills merge after the last image: {:.2f} s", merged->seconds);
+4 -37
View File
@@ -4,7 +4,6 @@
#pragma once
#include <array>
#include <atomic>
#include <condition_variable>
#include <exception>
#include <functional>
@@ -13,7 +12,6 @@
#include <map>
#include <memory>
#include <mutex>
#include <shared_mutex>
#include <vector>
#include "StillsPrePass.h"
@@ -26,7 +24,6 @@
#include "../common/ThreadSafeFIFO.h"
#include "../image_analysis/IndexAndRefine.h"
#include "../rugnux/StillsLiveMerge.h"
#include "../rugnux/StillsRunningGeometry.h"
// Serial frames analysed on the GPU as the DECTRIS ones are (MXAnalysisWithoutFPGA), after a pre-pass on
// the first PREPASS_IMAGES (RunStillsPrePass) has decided the beam centre and the mask they are analysed
@@ -44,17 +41,10 @@
// before Finish returns (GetMerged). The cell is the run's consensus at the first merge, as rugnux merges in
// it (the collection's cell picks the frames merged); the space group the collection's, else P1.
// merge = false leaves it out.
// Where the pre-pass refined the centre and distance on fewer frames than StillsRunningGeometry bundles, or
// could not refine them (too few of the frames it had indexed: a small image buffer, a serial hit rate, a chip
// whose first frames are empty), and a cell is given, they are refined again while the frames are analysed
// (StillsRunningGeometry, from the spots and lattices of the frames analysed, no image held), and the analysis
// switches to that geometry as each of GEOMETRY_ROUNDS fits is made. The frames analysed before the last
// switch are not merged (their reflections were predicted at an earlier geometry; the images are gone); where
// no fit can be made they are merged as analysed.
// Moves the analysis to a refined geometry during the run. The q range stays where it was, also one left to
// reach the detector's edge: the start message and the run's profile are binned on it.
void MoveStillsGeometry(DiffractionExperiment &experiment, float beam_x_pxl, float beam_y_pxl, float distance_mm);
// The geometry is fixed when the pre-pass ends: every frame is analysed and merged at it. Where the pre-pass could
// not refine it (too few of the frames it had indexed: a small image buffer, a serial hit rate, a chip whose first
// frames are empty; or no cell given), the header's distance stays for the whole run, with a warning in the log
// and the status.
class StillsGPUAnalysis {
public:
// Enough frames for the distance to be refined at a serial hit rate: of 100 frames of a collection indexed
@@ -147,29 +137,6 @@ private:
std::optional<StillsMergeOutput> merged;
std::future<void> merge_future;
// The geometry refined while the frames are analysed (see above), in GEOMETRY_ROUNDS rounds, each on the
// frames analysed at the last one's geometry (the spots of a frame are indexed at the geometry it was
// analysed at, and a start 2.6% off left the first fit 0.6% off). The analysis takes geometry_mutex shared
// for each frame; a switch takes it exclusive and raises geometry_epoch, which the analysis threads see to
// rebuild their analysis. geometry_state is PENDING until the last round, then SWITCHED or KEPT.
// frame_epoch: per frame number, 1 + the epoch it was analysed at while PENDING (0: after); the merge
// leaves out those analysed before the last switch.
constexpr static int GEOMETRY_ROUNDS = 2;
enum GeometryState : int { GEOMETRY_PENDING = 0, GEOMETRY_SWITCHED = 1, GEOMETRY_KEPT = 2 };
std::unique_ptr<StillsRunningGeometry> running_geometry;
int geometry_round = 0;
std::shared_mutex geometry_mutex;
// A switch waiting for geometry_mutex closes this gate to new frames, or the frames that keep taking the
// mutex shared would hold it off to the end of the collection.
std::mutex gate_mutex;
std::condition_variable gate_cv;
bool switch_waiting = false;
std::atomic<int> geometry_epoch{0};
std::atomic<int> geometry_state{GEOMETRY_KEPT};
std::vector<uint8_t> frame_epoch;
std::vector<std::unique_ptr<AzimuthalIntegrationMapping>> retired_mappings; // the analyses before the switch read them
void SwitchGeometry();
void StartPrePass();
void PrePassThread();
void AnalysisThread();
-2
View File
@@ -25,8 +25,6 @@ ADD_LIBRARY(JFJochRugnux STATIC
RugnuxLive.h
StillsLiveMerge.cpp
StillsLiveMerge.h
StillsRunningGeometry.cpp
StillsRunningGeometry.h
RugnuxRefine.cpp
RugnuxRefine.h
Rugnux.h
-95
View File
@@ -1,95 +0,0 @@
// SPDX-FileCopyrightText: 2026 Filip Leonarski, Paul Scherrer Institute <filip.leonarski@psi.ch>
// SPDX-License-Identifier: GPL-3.0-only
#include "StillsRunningGeometry.h"
#include <algorithm>
#include "../image_analysis/scale_merge/CellPopulation.h"
namespace {
// A frame the geometry refinement can use, as rugnux takes it from its sample: the spots on lattice 0 with
// their index, the others without.
GeomRefineFrame RefineFrame(const DataMessage &msg) {
GeomRefineFrame f;
f.lattice = *msg.indexing_lattice;
f.ordinal = static_cast<int32_t>(msg.number);
for (const auto &s: msg.spots) {
if (s.indexed && s.lattice == 0)
f.spots.push_back(GeomRefineSpot{s.x, s.y, static_cast<int32_t>(s.h), static_cast<int32_t>(s.k),
static_cast<int32_t>(s.l), s.panel});
else
f.unindexed.push_back(GeomRefineSpot{s.x, s.y, 0, 0, 0, s.panel});
}
return f;
}
}
StillsRunningGeometry::StillsRunningGeometry(const DiffractionExperiment &in_experiment, int in_nthreads)
: experiment(in_experiment), nthreads(in_nthreads) {}
bool StillsRunningGeometry::Add(const DataMessage &msg) {
std::optional<GeomRefineFrame> f;
if (msg.indexing_result.value_or(false) && msg.indexing_lattice.has_value())
f = RefineFrame(msg);
std::unique_lock ul(m);
if (done)
return false;
++offered;
if (f && f->spots.size() >= 6) {
++indexed;
strong += static_cast<int>(f->spots.size()) >= MIN_STRONG_SPOTS;
frames.push_back(std::move(*f));
}
done = static_cast<int>(strong) >= STRONG_FRAMES || indexed >= INDEXED_BUDGET;
return done;
}
int64_t StillsRunningGeometry::Offered() const {
std::unique_lock ul(m);
return offered;
}
size_t StillsRunningGeometry::Strong() const {
std::unique_lock ul(m);
return strong;
}
std::optional<GeometryRefinerResult> StillsRunningGeometry::Refine() {
std::unique_lock ul(m);
const auto cell = experiment.GetUnitCell();
if (!cell)
return {};
// As RugnuxStills.cpp's sample: only the frames of the given cell's population are fitted, the population
// told from the strong frames.
const double given = CrystalLattice(*cell).CalcVolume();
std::vector<UnitCell> strong_cells;
for (const auto &f: frames)
if (static_cast<int>(f.spots.size()) >= MIN_STRONG_SPOTS)
strong_cells.push_back(f.lattice.GetUnitCell());
if (strong_cells.empty())
return {};
const auto [low, high] = GivenCellPopulation(strong_cells, *cell);
std::vector<GeomRefineFrame> bundle;
for (auto &f: frames) {
const double r = f.lattice.CalcVolume() / given - 1.0;
if (static_cast<int>(f.spots.size()) >= MIN_STRONG_SPOTS && r >= low && r <= high)
bundle.push_back(std::move(f));
}
frames.clear();
std::sort(bundle.begin(), bundle.end(), [](const GeomRefineFrame &a, const GeomRefineFrame &b) {
if (a.spots.size() != b.spots.size()) return a.spots.size() > b.spots.size();
return a.ordinal < b.ordinal;
});
if (static_cast<int>(bundle.size()) > STRONG_FRAMES)
bundle.resize(STRONG_FRAMES);
GeometryRefinerSettings settings;
settings.crystal_system = experiment.GetCrystalSystem();
settings.num_threads = nthreads;
settings.cell_reg_coeff = 0.0; // cell and distance both free, as rugnux fits them
const GeometryRefinerResult r = RefineGlobalGeometry(experiment.GetDiffractionGeometry(), *cell, bundle, settings);
if (!r.ok)
return {};
return r;
}
-47
View File
@@ -1,47 +0,0 @@
// SPDX-FileCopyrightText: 2026 Filip Leonarski, Paul Scherrer Institute <filip.leonarski@psi.ch>
// SPDX-License-Identifier: GPL-3.0-only
#pragma once
#include <mutex>
#include <optional>
#include <vector>
#include "../common/DiffractionExperiment.h"
#include "../common/JFJochMessages.h"
#include "../image_analysis/geom_refinement/GeometryRefiner.h"
// The detector centre and distance of a serial collection refined from frames as they are analysed, holding
// no image: the live path's geometry where its pre-pass could not refine it (too few of the frames the image
// buffer holds indexed). The frames' indexed spots and lattices are kept until STRONG_FRAMES strong frames
// (at least MIN_STRONG_SPOTS indexed spots) are in, or INDEXED_BUDGET frames indexed (empty frames - a chip's
// wells without a crystal - do not count); then the strongest
// are bundled and the geometry refined once, as rugnux refines it on its sample (cell and distance free, the
// detector untilted), on the frames of the given cell's population. A fixed cost: one bundle, one fit.
// Thread-safe.
class StillsRunningGeometry {
public:
constexpr static int STRONG_FRAMES = 200;
constexpr static int MIN_STRONG_SPOTS = 20;
constexpr static int64_t INDEXED_BUDGET = 5000;
// `experiment`: the geometry the frames are analysed at, the cell given.
explicit StillsRunningGeometry(const DiffractionExperiment &experiment, int nthreads);
// A frame as analysed. True once: for the frame after which Refine is to be called.
bool Add(const DataMessage &msg);
// The refined geometry, where the fit ran.
std::optional<GeometryRefinerResult> Refine();
// Frames offered and frames kept, for the log
[[nodiscard]] int64_t Offered() const;
[[nodiscard]] size_t Strong() const;
private:
const DiffractionExperiment experiment;
const int nthreads;
mutable std::mutex m;
std::vector<GeomRefineFrame> frames;
int64_t offered = 0, indexed = 0;
size_t strong = 0;
bool done = false;
};
-27
View File
@@ -125,30 +125,3 @@ TEST_CASE("StillsGPUAnalysis_FailingFrameDoesNotStopAnalysis", "[StillsGPUAnalys
CHECK(s.steps.back().duration_s);
CHECK(buffer.GetStatus().preparation_slots == 0);
}
TEST_CASE("StillsGPUAnalysis_MovedGeometryKeepsTheProfileBins", "[StillsGPUAnalysis]") {
DiffractionExperiment x = StillsExperiment();
const PixelMask mask(x);
const AzimuthalIntegrationMapping before(x, mask);
AzimuthalIntegrationProfile run_profile(before);
// Left to itself the high q follows the detector's edge, and with it the number of bins.
DiffractionExperiment unfrozen = x;
unfrozen.DetectorDistance_mm(120);
REQUIRE(AzimuthalIntegrationMapping(unfrozen, mask).GetBinNumber() != before.GetBinNumber());
MoveStillsGeometry(x, 510, 260, 120);
CHECK(x.GetDetectorDistance_mm() == Catch::Approx(120));
CHECK(x.GetBeamX_pxl() == Catch::Approx(510));
CHECK(x.GetBeamY_pxl() == Catch::Approx(260));
const AzimuthalIntegrationMapping after(x, mask);
CHECK(after.GetBinNumber() == before.GetBinNumber());
CHECK(after.GetBinToQ() == before.GetBinToQ());
const AzimuthalIntegrationProfile frame_profile(after);
CHECK_NOTHROW(run_profile += frame_profile);
// A second move stays on the same bins.
MoveStillsGeometry(x, 505, 255, 90);
CHECK(AzimuthalIntegrationMapping(x, mask).GetBinToQ() == before.GetBinToQ());
}