Making the worker streams non-blocking removed the implicit ordering that the
constructors were still relying on. Each engine uploads its static inputs - the
pixel mask, the pixel-to-bin map, the corrections, the ROI map - with a blocking
NULL-stream cudaMemcpy, and then reads them from kernels on its own stream. A
pageable host-to-device cudaMemcpy returns once the source has been staged, with
the DMA still in flight, and a non-blocking stream no longer waits for the NULL
stream. The failure mode is a silently unapplied mask or a stale mapping, not a
crash, so it would not have announced itself.
Put them on the stream the engine already owns, and synchronise once at the end of
the constructor - that is required for the preprocessor, whose source is a local
vector, and leaves the others settled rather than in flight for the cost of one
one-time sync. The GPU spot-finder test uploaded its image the same way.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Leaving it at G = 1 looked like the conservative choice and is the more damaging
of the two errors. The per-image scale enters as rlp/(partiality*G) and multiplies
intensity and sigma alike, so substituting 1 for a scale that was really 1/200 of
the run median puts the intensities in 200x too low with sigmas 200x too low too -
1/G^2 times the weight they deserve. The merge cannot defend itself against that,
because the number that is wrong is the number the weight is built from. And if
the collapsed value was instead a failed fit, G = 1 merges the image mis-scaled by
an unknown factor. Per-crystal scales on serial stills genuinely span orders of
magnitude, unlike frames of one rotation sweep, so both readings are live.
An image whose scale is not believable has no usable scale. Write NaN into its
image_scale_corr, which every merge path already skips on, so it drops out of the
merged intensities, the error model and the statistics consistently.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The "keep what the crystal came in with" gate required std::isfinite(cc) before it
would reject, so a refined model whose CC could not be measured at all was adopted.
ImageReferenceCC returns NaN when fewer than 20 reflections clear the partiality
cut - which is exactly what a refinement that collapsed the partialities produces,
since the cut is on the partialities it just rewrote. The gate therefore failed
open on precisely the crystals it exists to catch, and wrote the NaN into
image_scale_cc, on which --min-image-cc then drops the image from the merge, the
error model and the statistics.
Treat a CC that cannot be measured as worse than one that can, so the crystal is
put back exactly as it arrived.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
"Analyze dataset" cleared the stored cell and space group unless "Use the stored
unit cell / space group" was ticked, and that checkbox defaulted off. But the
settings panel writes the user's own cell and space group onto the experiment, so
a cell typed into the panel was discarded too - while the checkbox label said
"stored", implying it came from the file.
It also contradicted the dialog next to it: "Refine geometry (stills)" is offered
and default-ticked precisely because a cell is present, and the run then removed
that cell. The default dialog state on a stills dataset with a known cell ran the
bundle adjustment with nothing to anchor on and dropped indexing off ffbidx, which
needs a cell, onto de-novo FFT.
Drop the checkbox and take the crystal from the panel, which already has exactly
the right semantics: "Unit cell known" ticked writes the cell and group, unticked
clears both, and a space group of 0 means none. So ticked = -C/-S, unticked =
bare rugnux, and what a run will use is always what is on screen. That also keeps
the copied command line honest, since RugnuxCommandLine emits -C/-S from the same
experiment. The panel is refilled from the file when one is opened, so a finished
job's _process.h5 becoming the active snapshot now shows its group and can be
cleared, instead of silently pinning every later run to it.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
--polarization was applied with the other geometry overrides, but
configure_offline_output runs afterwards and calls ApplyRugnuxExperimentDefaults,
which sets the polarization factor unconditionally. Every full-analysis run used
0.99 whatever was asked for, so the Lp correction was wrong at a beamline with
different polarization. Apply it after the defaults instead, and stop claiming in
RugnuxDefaults.h that nothing here is user-selectable.
--scale built a bare ScalingSettings and re-derived the rotation/stills split by
hand rather than calling RugnuxDefaultScalingSettings, which is what the split was
factored out for. It got scale-fulls, smooth-G, min-captured-fraction and outlier
rejection right and dropped CaptureUncertaintyCoeff on the floor: 1.0 in the
pipeline, 0.0 here. So re-scaling a rotation _process.h5 gave different sigmas and
ISa than the run that wrote it - the exact failure the block's own comment says it
exists to prevent. Start from the shared defaults and apply the overrides on top,
which also picks up --mosaicity and --search-min-zeta, and let -C bind here too.
REJECT_OUTLIERS_DEFAULT_NSIGMA had no reader left afterwards.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The kernel guards against 2*max_hkl+1 and maps thread i to h = i - max_hkl, but
the host launched a grid sized 2*max_hkl. The h = k = l = +max_hkl planes were
therefore never launched while -max_hkl was, so the GPU predicted an asymmetric
subset of what the CPU loop (inclusive on both ends) does. The same bug was fixed
on the stills twin when the whole hkl range moved to the GPU; the rotation
predictor kept the old expression.
It only bites where the cell actually reaches |h| = 100 inside d_min - a ~150 A
axis at 1.5 A - so most data never noticed. Over the 33-crystal rotation battery
29 crystals are bit-identical and 4 gain observations, all of them large-cell or
high-resolution: +8519, +4693, +901 and +758 observations, with the high-shell
CC1/2 up 15.0->15.1%, 52.0->52.2%, 76.3->76.6% and 51.6->52.1%. Nothing is lost
anywhere, and R-meas and ISa move by at most 0.01.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The first pass sampled min(n, max(refine_frames * 50, 8000)) images to keep the
200 strongest, so any serial run of 10000 frames or fewer indexed every frame
TWICE - and 99.3% of the pass was that sampling, the bundle adjust itself taking
0.24 s. The budget is sized for a low-hit-rate dataset; on data that indexes well
almost all of it was wasted.
Stop once four times the bundle size has been found, which still leaves the
"strongest N" selection a real pool and still spans the run, because the sample
is equally spaced. On a lysozyme jet dataset that is 1500 frames examined instead
of 4000, the same 200 bundled, and the same refined geometry - beam and distance
to the pixel, cell to 0.01 A. Warm cache: the pass drops 33 s -> 9.1 s and the
whole run 65 s -> 35 s, with CC1/2 and R-meas unchanged inside replicate noise.
Stills only - rotation has its own two-pass and returns from this function early.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Every per-thread stream was created with cudaStreamDefault, and the 20 MB raw
image upload went to the legacy NULL stream. A NULL-stream operation implicitly
synchronises with every blocking stream in the process, so with one engine per
worker thread no two workers' GPU work could ever overlap - the whole GPU
pipeline ran serially however many threads were asked for.
Create the streams non-blocking and put the upload on the engine's own stream.
Measured on 2000 serial stills, interleaved, medians of three: 24.6 -> 19.2 s at
-N 32 (-22%), 32.7 -> 21.0 s at -N 16 (-36%), CPU utilisation 436-570% -> 723-859%.
Output bit-identical - same observations, uniques, completeness, R-meas, CC1/2,
error model and cell. The stream is synchronised at the end of the same function,
so the ordering the code relies on is unchanged.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
ScaleOnTheFly's collapsed-scale guard ran, and then StillsPartialityRefine
re-fitted every crystal's scale with no floor and adopted it unconditionally
whenever the image had no prior CC - which is exactly the state the guard leaves
behind. So the guard was protecting almost nothing. Measured on a lysozyme jet
dataset: of 367 images it left unscaled, only 8 were still unscaled in the
output, and 53 reached the merge at or below a fiftieth of the run median, the
worst at a 4525th; on a second run of the same sample, 297 images, worst at a
75000th. Those intensities are what the merge saw - up to 94x too high in the
written file.
Run the guard again on the refined scales. It now reports 367 then 51 on that
dataset, and the merge improves: R-meas 117.4 -> 110.8%, CC1/2 96.3 -> 96.5%.
Measurements confirm the rest of the guard is right as it stands: 0.02 is ~5x
below the lowest scale ever seen on an image that correlates with the merge
(no image at CC >= 0.4 falls below a tenth of the median), and leaving the image
at G = 1 beats both dropping it and replacing its scale with the median - on a
run where 30% of images are affected, dropping costs 1.8 CC1/2 and 29%
multiplicity.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Raising the per-image reflection limit to 65536 for offline reprocessing also
raised the image-buffer headroom derived from it, and that headroom divides a
FIXED total buffer - so every slot grew from compressed+4 MB to compressed+16.7
MB and the receiver's slot count, i.e. how much of a burst it can absorb, fell by
about three. Online never needed the raised limit: measured on three serial
stills datasets the worst frame predicts 1380 reflections, 14% of even the old
cap.
So split them, the same way the geometry refinement's stopping rule is split:
online keeps the transport-sized 10000, offline gets the full 65536, and the
buffer headroom derives from the online one. Both still come from BraggPrediction
so the cap, the prediction and the headroom cannot drift apart.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
RefineOne re-measured the image's correlation to the reference after writing the
refined partialities - because --min-image-cc drops images by it - and then
ignored what it measured. A crystal the tilt model suits worse than the fixed
partiality it replaces kept the refined model anyway, and the refinement is on by
default. Compare against the CC the crystal arrived with and put it back
untouched when the refinement does not improve it, which is the same state a
crystal with too few reflections to fit ends in.
Also four things noted in review and left until now: AdaptiveThresholdTest.cpp
was listed twice in the test target, AdaptiveThreshold.h was the one header in
image_analysis/spot_finding not in its library's source list, CLAUDE.md said
update_version.sh rewrites VERSION when it only reads it, and the CHANGELOG did
not mention that image_scale_b is gone from the plot_type enum - which breaks a
client that asks for that plot.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The per-image refinement stopped on a wall-clock budget (40 ms, and 20 ms for
the rotation-only extra pass). Online that is exactly right - the budget is real
and an image that overruns it costs the acquisition. Offline it means the same
file refines to a different lattice depending on what else the machine was doing
at the time, which is not a property reprocessing should have.
Bound it by iteration count instead when the caller is offline. IndexAndRefine
takes the workflow as a constructor argument: the receiver asks for the
wall-clock bound, rugnux and the viewer get the reproducible one. 50 iterations
is Ceres' own default; the per-image problem converges well inside it, so it
bounds the pathological case rather than the normal one - measured on five
battery crystals, every number is unchanged from the timed version.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The second pass re-indexes de novo and can land in a different setting from the
first - most often on the PRIMITIVE sub-cell of a centred lattice. The reindex
that exists to undo that declines when the metric does not match, and the code
then went on to stamp pass 1's group onto the cell regardless.
That is not a small error. A C-centred group on an already-primitive cell means
the centring absence rule removes half the reflections that genuinely exist, so
the merge holds more unique reflections than its own cell can - measured here as
"117% complete" with CC1/2 0.62, against the first pass's 92.6% and 0.98, on a
cell of exactly half the C-centred volume.
Detect the conflict where it happens and feed it to the pass-2 credibility guard
rather than acting on it locally: letting the pass re-search its own group
instead produced a P1 answer on a crystal XDS and the first pass both call C2,
which is a worse outcome than simply not trusting the pass. The completeness and
CC1/2 tests stay as the symptom-side net.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Digging into the selection logic showed the caps were not deciding the science -
the two-pass geometry post-refinement was, and the caps only fed it randomness.
Caps. The prediction buffer now grows to whatever a frame predicts instead of
keeping an arbitrary subset of it, and the per-image reflection limit is raised
to 65536, with the image-buffer transport headroom derived from the same
constant so the two cannot drift. Measured: bit-identical output on five battery
crystals, because a normal cell never approached the old limits - only a large
cell (~2.8e6 A^3, ~30000-44000 predictions per frame) ever did.
Pass-2 guard. The refined pass is normally the better answer, which is why it is
the canonical output, but it was adopted whatever it produced. On that same
crystal it merged more unique reflections than its own cell can hold -
completeness "117%", which is arithmetically impossible - while the header-
geometry pass sat at 92.6% and CC1/2 0.98. Compare the two and, when the refined
pass is not credible, go back to the header geometry and re-run so the canonical
files are the ones that are kept. Both bounds are set where only a failure
reaches them.
Together on that crystal: 111639 unique against XDS's 118730 (was 88000-99000
and different every run), CC1/2 98.0% (was 96.9-97.7%), ISa 8.54, and two runs
now agree bit for bit.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Found while chasing a 12% run-to-run spread in the merged reflection count of one
crystal. The GPU kernels claim output slots with an atomicAdd and, on overflow,
undid the increment with an atomicSub - so the counter saturated at the capacity
and the host could not tell a full buffer from an overflowing one. Which
reflections survived was then decided by CUDA block scheduling and changed every
run. Measured on that dataset: every frame predicts 23000-44000 against a 20000
buffer, and the spread reached the merged output (161591 / 165193 / 166110 /
166479 unique across four runs of the same command). Single-threaded runs diverge
too - this is entirely GPU-side.
Stop clamping the counter, so the true number predicted reaches the host, and
warn once per predictor when it exceeds the buffer. Which reflections are kept is
unchanged: making that reproducible means deciding what to keep when a frame
predicts more than the pipeline carries, and the obvious answers are worse - the
capacity is not the real limit, kPredictionOutput (10000, selected by smallest
excitation error) is, and on this crystal both a bigger buffer and a strided
selection collapse the merge, because the rotation combine rebuilds fulls from
exactly the partials that a smallest-excitation-error cut throws away.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
"Analyze dataset" and `rugnux` with no options are two front ends onto the same
library and are meant to agree, but they decided their defaults separately and
the two lists had drifted. The viewer was missing:
- the de-novo starting point. The CLI discards the cell and space group stored
in the input file before it does anything; the viewer left them on the
experiment. Rugnux only searches for a space group when none is set, so the
search was skipped entirely and the stored group was reported straight back.
That is self-reinforcing: a finished job's own _process.h5 becomes the active
snapshot, so a run that ended in P1 pinned every later run to P1 - which is
what "lysozyme keeps coming out P1 in the viewer" was.
- the polarization factor, so the Lp correction was omitted altogether. The
missing factor is azimuthal and intensity-proportional, and symmetry mates sit
at the same 2-theta but different azimuth - the exact "unequal intensities
forced together" signature the space-group search vetoes as pseudo-symmetry,
which can land a genuinely de-novo run in P1 on its own.
- five rotation scaling defaults: the smooth-G range, the minimum captured
fraction, the capture-aware sigma, outlier rejection and --search-min-zeta.
The CLI's own comments tie the captured-fraction default to a crystal
recovering its true space group instead of P1.
Put the policy in one place (RugnuxDefaults) and have both front ends start from
it. The CLI now takes its defaults from there and applies user options on top;
its output is unchanged, verified bit-for-bit on four battery crystals.
The stored cell/group is still available: the job dialog offers "Use the stored
unit cell / space group", off by default, shown only when the file has one.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The error model says sigma -> b*I for strong reflections, so merged I/sigma
flattens off at 1/b - the number reported as ISa. Plotting I/sigma against I
with that asymptote drawn on it is what shows whether the reported ISa
describes the data or comes from a degenerate fit, which nothing in the window
could show before. A third page next to the per-shell plot and table.
The merge carries a few thousand strided (I, sigma) pairs to the viewer for it -
a shape, not a reflection list; the reflections themselves are in the .mtz/.cif.
The hero row gains the Wilson B and the radiation-damage Delta-B. Both were
already computed and already in MergeStatistics, so they only needed showing.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The protection against a per-frame scale collapsing toward zero lived inside
ComputeSmoothGWindow, so it only existed when smooth-G did: --smooth-g=0, a
dataset whose oscillation width is unknown, and any caller that never sets a
smoothing range - the viewer among them - merged with no guard at all. A
collapsed G multiplies that frame's intensities by 1/G and its sigmas by the
same factor, so nothing downstream can see it; the merge's n-sigma cut scales
with the number that is wrong.
Pull it out into ReplaceCollapsedScales, called unconditionally right after the
partial scaling loop, and let the smooth-G window assume what it now guarantees
instead of computing its own median and floor.
The fulls guard built its median from every frame including those never fitted -
those sit at the combine's corr = 1, so a run with many unfitted frames dragged
the median toward 1 and the floor with it. It also reported the absolute
amplification where the message says "below the run median".
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Also corrects the reader bit-depth entry: taking the depth from the file was
the wrong fix for the 32-bit EIGER2 case and broke 8-bit files.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Eight pages under docs/python_client/docs describe schemas that appear nowhere
in jfjoch_api.yaml and are linked from no index - left behind because the
regeneration step never cleared python-client/, which these are copied from.
That is fixed in update_version.sh; this removes what accumulated.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
None of this has a reader:
- ScalingSettings::scaling_regularize and its setter/getter
- ScaleOnTheFlyResult::succesful (never set) and ::time_s (set, never read),
with the timing that only fed the latter
- JFJochImage::last_fit_viewport_ (written twice, read nowhere) and the
comment claiming the retry uses it - the retry keys off initial_fit_done_
- JFJochDiffractionImage::ice_ring_width_Q_recipA, and a QtConcurrent include
in a file that uses none
- an unused gemmi::Op accumulator in the spindle-angle helper
- <random> in Merge.{h,cpp}, from before the half-set split became a hash
- an orphaned comment describing the Ceres B-factor residual deleted in
014e43a4c, and two trailing comments that had collided on one line
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The script now removes the generated C++ model, the frontend client and the
published python docs before regenerating, but not python-client/ itself - and
docs/python_client/docs is filled by copying that directory. So a schema dropped
from the API kept its generated model in the PyPI package and its .md page in
the published docs, linked from no index. Eight such pages are in the tree
today, JfjochSettingsSsl among them.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Adaptive spot detection became the default for rotation data as well in
6f4917dce; RUGNUX.md still said rotation kept the fixed-threshold finder, which
is also the opposite of what the usage message and CPU_DATA_ANALYSIS say.
--search-min-zeta had no entry in the option tables at all.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
centerAt returns early while the window is hidden, and nothing replays the last
position when it comes back, so re-opening the magnifier showed whatever region
the cursor was over when it was closed - with current pixels, which makes it
look like a live view of the wrong place. Remember the position while hidden and
apply it on show.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
drawPixelLabels took its range from the whole viewport rather than from the
exposed rect it was given, so a 200x40 px hover repaint still walked up to 5000
cells doing mapFromScene + QImage::pixel + drawText for each, only to have the
result clipped away. With hover feedback now rate-limited to 15 Hz that ran ~75
times a second at high zoom, against once per overlay rebuild before the
rendering rework. Intersect with the exposed rect; same in the magnifier.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The "d = ... A" readout moved from a scene item flagged
ItemIgnoresTransformations to a fixed viewport position painted in
drawForeground, which is what made a hover update dirty a small rect instead of
the whole viewport. But QGraphicsView pans by blitting the viewport: the painted
text is shifted along with the image and left there, and the pending update for
its old position is translated away too, so dragging the image smears ghost
copies of the readout across the corner. Dirty the old and the new rect when the
view scrolls - still a couple of hundred pixels, not the viewport.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
With merging on, the _process.h5 is skipped because the merged reflections are
the wanted output and that file is large (113 MB for 200 images here). But if
nothing indexes there are no merged reflections either, so the run finished
successfully having written no file at all - the one case where the user most
needs something to look at.
Write it in that case. The per-image messages have already gone past unwritten,
so this carries the dataset metadata, the mask, the azimuthal profile and the
summary scalars rather than the full per-image tables - and it is small for the
same reason it is needed (124 kB on a zero-index run). A run that does index is
untouched.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
main had no enclosing try/catch, but plenty of ordinary input reaches a setter
that throws: --polarization 2, --detector-distance 0, -q 0, --azim-max-q 20,
--scale combined with a reference MTZ, and every failure inside the pipeline
itself. All of them ended as "terminate called after throwing an instance of
'JFJochException'" and exit 134, with the message nowhere to be seen. Move the
body into RunRugnux and let main report what was thrown, exit 1.
Two options also still bypassed the numeric parser that exists to prevent this:
--scaling-high-resolution used atof, which turns a typo into 0 and then throws
from the setter, and --integration-radius used std::stof, which throws on
non-numeric input.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The Ceres summary was discarded, so a solve that diverged or aborted left its
last iterate in psi and that tilt was written onto the partiality of every
reflection of the crystal. Restore the tilt the crystal came in with and stop
refining it; the scale fit alone is still a usable model, which is what the
other three early returns in this function fall back to.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
SolveScaleIRLS returns whatever it converged to and both writers accept any
G > 0, so a fit that collapsed to ~1e-3 multiplies that image's intensities by
a thousand. Nothing downstream notices, because the sigmas are multiplied by the
same factor and the merge's n-sigma outlier test is therefore blind to it - only
a total collapse self-heals, by overflowing corr to inf.
The rotation path refuses a per-frame scale this far below its neighbours; the
stills path had no guard. Judge each image against the median of the images that
did scale, and put a collapsed one back to G = 1 - the same state as an image
with too few reflections to fit.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Two faults in the same block, both of which let a search pass corrupt the
production merge that follows it.
The device's corr was only copied back to the host for the diagnostic dump, but
the |zeta| filter runs on the host and then uploaded the whole host array - so on
a CUDA build it wrote the values ingested BEFORE scaling over the scaled and
smoothed corr the device had just computed. With the rotation default
--search-min-zeta 0.85 that means the space-group search was deciding the
symmetry from an unscaled merge. Copy corr back first, and upload once after
both filters instead.
Zeroing corr also has no owner: it is how an observation leaves the merge, but
the only thing that ever rewrites it is the scaling loop, which skips frames it
cannot fit. A frame left with too few well-measured reflections therefore kept
its dropped observations at zero for the rest of the object's life - and the
final production merge re-uses the same object without re-ingesting. Snapshot
corr before the filters and restore it at the start of the next pass, so each
pass decides for itself and the final merge keeps everything, as documented.
The frame rejection (--min-image-cc) is now applied on the host for both paths;
its separate device path did nothing whenever the CPU combine was in use.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The existing cases plant blobs at 200 on a background of 8..12, so any threshold
between 12 and 200 passes them - replacing RingThreshold with a constant leaves
them all green. Two cases that do not:
- the CPU threshold has to track the background: a frame and the same frame
scaled ten times must give the same spots, with a pixel a few sigma above the
background staying unfound in both. A constant threshold, or one that drops
the sigma term, fails one scale or the other.
- the GPU engine has to agree with itself across runs, which is what the ring
sums being order-independent buys.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The ring sigma is the cancelling difference sum2/n - m^2, and both sums were
float accumulated by atomics whose order is arbitrary. Two costs: the
cancellation left only ~4 digits in the variance, and the ordering moved the
resulting threshold by ~0.05 counts between runs - enough to flip a pixel
sitting on the hard "value >= threshold" test, and with it a connected
component's size. So the GPU engine did not reproduce the CPU one and did not
reproduce itself.
Only the accumulators that span blocks are widened. The per-block staging stays
float, because a block contributes a few dozen similar-magnitude pixels to a
ring and there is nothing to lose there - that also keeps the shared-memory
footprint of the hot loop, and hence its occupancy, exactly as it was: measured
on a 4.5 MP frame, 0.960 vs 0.966 ms/frame (40.9x over the CPU path, unchanged).
finalize_rings now does the cancellation in double and rounds to float last,
which is what AdaptiveSpotFinderCPU::AccumulateRings does.
The device properties are also read from the current device rather than device
0; callers round-robin engines across GPUs.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
b81c6f00b took the container depth from bit_depth_image in the file, which fixed
32-bit EIGER2 files but got the general case wrong: the reader converts every
image to SIGNED int32 (PixelSigned(true) a few lines up), while bit_depth_image
describes an unsigned container, and GetOverflow() combines the two. So a 16-bit
file still capped at INT16_MAX rather than 65535, and an 8-bit file newly capped
at 127 - flagging counts 127..254 as saturated, which drops the whole reflection
at the integration accept gate.
Declare 32 bits, matching what the reader actually hands out. The saturation cap
then comes from the file's own saturation_value, which is what it is for, and
the error value reported for a read dataset becomes INT32_MIN - the sentinel the
reader really uses. A file whose bit_depth_image is not 8/16/32 also stops
throwing on open.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The host sized the grid from 2*max_hkl while the kernel guards against
2*max_hkl+1, so whenever the rounded-up grid landed exactly on 2*max_hkl threads
(max_hkl a multiple of 4, with the 8x8x8 block) the h = +max_hkl plane was never
launched. The CPU loop runs -max_hkl..+max_hkl inclusive, so the GPU predicted a
strict subset.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
calc_std uses the cancellation-prone (sum2 - sum^2/n) form on float accumulators
summed over millions of pixels, so a flat ring - true variance near zero - comes
out negative as often as positive and GetStd() returns NaN. Both adaptive
spot-finder ring accumulators already floor this at zero; this one did not, and
1a0774eed made sum2 correct, so the path is now actually exercised.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Making the high limit optional removed the implicit upper bound on the low one
(it used to follow from high <= maxQ and high > low), so --azim-min-q 50 with no
maximum is accepted and ResolveHighQ then calls std::clamp with its lower bound
above its upper bound, which is undefined.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Every other reader of spot_finding.high_resolution_limit spells "unset" as
value_or(0) and compares, so 0 and nullopt are interchangeable - except in
SpotAnalyze, which passed the 0 straight to ResolutionShells and threw
"Resolution must be above zero" on every image. Reachable over the REST API,
where 0 is the natural way to say "no limit" and the settings check lets it
through; the rugnux CLI already maps 0 to unset before this point.
While here, check that a limit that IS set is finite regardless of its sign -
NaN fails the > 0 test and was skipping validation entirely.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The two-pass rotation run reinstates pass-1's space group for pass 2's merge, so
pass 2 skips the search block entirely - and with it the flag that records that
the Laue class was CHOSEN by the search rather than given. The canonical output
therefore printed the plain "no twinning: the Laue class is holohedral, so no
merohedral twin law exists", which is exactly the circular conclusion the flag
was added to replace; only the throwaway _01 output carried the caveat.
The text is written per pass, inside RunPipeline, so the flag has to travel with
prepass_merge_sg_ rather than being patched onto the returned result.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The two-arm search is meant to use the Lorentz-filtered merge for the point
group only - systematic absences live in the weak reflections a |zeta| cut
removes, and reading them off the filtered arm is what cost four crystals their
screw axes. That is what the code comment and CPU_DATA_ANALYSIS both say, but
the filtered-arm-wins branch kept its whole result, screws and centering
included.
Let a search be pinned to a point group decided elsewhere (fixed_point_group)
and re-run Stage B on the all-observation merge when the filtered arm rescues
the point group. The point group is passed as its symmorphic representative, not
by name: gemmi calls both P321 and P312 "32". Reporting that representative also
lets the ambiguity check see two arms that disagree about which 2-folds are real
- by name they looked identical - and the advice it prints now names a space
group -S can actually be given.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The row-relative cut scales the "too strong to be absent" threshold by the axial
row's own median E^2, floored at the plain value - so it can only raise the bar,
and a row whose control class holds a single strong reflection sets it from that
one reflection. That direction invents screws: a genuine 4_2 whose 00l happen to
be observed only at l=4n reads its l=4n+2 reflections as absent and ranks
4_1/4_3 above the truth. Require three controls before the row may set the
scale; below that the row keeps the plain cut.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
A candidate can have several confirmed subgroups of the same order - 422 has
both 4 and 222 - and on a twinned crystal the rival is not a harmless
alternative: a P4 crystal twinned by 2[100] has 222 confirmed too, and 222
CONTAINS the twin laws, so its own merge b is already ballooned. The H test
already answers to every tied parent; the systematic-b veto and rescue took
whichever one the enumeration happened to list first (222 before 4, by space-
group number), which disabled the veto on exactly the case it exists for. Take
the smallest parent b, which is the conservative direction for both tests.
The refusal message also quoted the raw parent b rather than the floored value
the veto actually compared against.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The cap that refuses to report an impossible ISa was zeroing the asymptotic b
itself, and that same value is the floor passed to SigmaWithSystematicFloor -
where zero means "no floor". So on the degenerate low-multiplicity fit the guard
is written for, instead of capping merged I/sigma at 100 it removed the cap
entirely. Report the asymptote as unmeasured, keep the fitted value for the
floor.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
6be94f2be stopped subtracting a negative intensity's Poisson term from the
background variance, but only in the host Combine(). The CUDA combine is the
path that actually runs: Run() selects it whenever a device is present and no
observation dump was asked for, so the correction never took effect on a normal
run, and a --dump-observations run merged differently from a normal one - the
two are meant to be identical.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The document had accumulated development narrative: what was tried and
rejected, which datasets a change rescued or cost, measured percentages
from test batteries. That belongs in commit messages and reports, not in
a reference description of the pipeline - it dates quickly, and a reader
looking up what an algorithm does has to sort it out from how it came to
be.
Removed throughout, keeping the algorithmic content and the design
reasoning that explains a choice on principle:
* 3.2 the whole paragraph justifying the rotation spot-finding limit
from battery measurements, and the CPU-vs-GPU per-frame timings;
* 3.3 "a significance/z-score was considered but is uninformative";
* 7.4 / 7.5 the comparisons to a robust loss and to joint refinement
as approaches that had failed;
* 9.2 the R_meas / CC1/2 outcomes attributed to the trimmed-mean
background;
* 9.3 "per-detector-region and crystal-anisotropy profiles were
evaluated and add nothing";
* 10.2 the stills tilt "succeeds where a freely-fitted width
collapses";
* 10.5 the CC_anom argument, trimmed to why the statistic behaves as
it does;
* 10.6 the survey of per-frame correlation medians across datasets;
* 13 the space-group bullet, restructured into the three gates it
actually applies, dropping the dataset anecdotes;
* 14.2 the free-form per-shell rescale, stated as a design choice
rather than an experiment.
Section 13's space-group text was one 20-line paragraph; it is now a
numbered list of the three tests, which is what the code does.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Rotation kept a 1.5 A high-resolution limit for spot finding on the
strength of one indexing-rate measurement (100.0 -> 96.8% on a strong
set). Measured properly, over the whole 33-crystal rotation battery, it
does not earn its place:
* no space-group decision changes - the same 30/33 agree with XDS, and
the three that differ are the known pre-existing cases, unchanged;
* 29 of 33 crystals are identical to the digit - same indexing rate,
R_meas, CC1/2, ISa. The limit was doing nothing on the large
majority;
* where it does bite, the limit is the worse setting. The one crystal
that loses appreciable indexing rate without it (99.50 -> 94.22%)
comes back with lower R_meas (29.4 -> 28.1), higher high-resolution
CC1/2 (27.9 -> 29.1) and higher ISa (5.77 -> 6.17). Another loses
0.4% of frames and gains 2.8 points of CC1/2_hi. Fewer frames
indexed, better data from them;
* runtime is unchanged (16m46s vs 17m32s over the battery).
So the indexing-rate cost is real but does not carry through to the
merged data, which is what the limit was protecting. Unset now means "as
far as the detector reaches" for rotation as well as stills;
--spot-high-resolution still sets a limit for weak, high-background data
where the extra high-resolution spots are genuinely noise.
This also removes the flag that distinguished "the user asked for no
limit" from "the user said nothing" - with no rotation default left,
both mean the same thing. While rewriting the comment block, corrects
its neighbouring claim that rotation keeps the fixed-threshold finder;
adaptive detection has been the default for both workflows since
6f4917dce.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Ten places where the document and the implementation had drifted apart.
Each was checked against the source before rewriting:
* 7.5 rotation post-refinement: it is TWO separate cross-validated
steps (cell+axis from the angles, then distance+beam from the
positions with the cell fixed), not one joint fit against the merged
fulls; the held-out split is an hkl hash, not a frame split; the
bounds are +-5% on distance and +-15 px on the beam, not "under
~1%"; and only the distance and beam centre reach pass 2, which
re-indexes de novo.
* 9.2 the trimmed-mean background: it is computed in the shared
background pass, so it DOES apply to --integrator boxsum. Only the
broadband sigma-clip is excluded. The section previously said both,
contradicting itself two paragraphs apart.
* 9.3 per-reflection profile rebuild, sub-pixel centring and radial
elongation are gaussian-only; the empirical profile keeps the fixed
per-shell grid and is accumulated on rounded predicted positions,
not centroids.
* 10.5 the asymptotic ISa and the b_ISa sigma floor are rotation-only;
stills report 1/b and floor with the whole-range b.
* 8.4 centering absences are applied only when the user fixes the
space group - de novo, prediction runs in P so the search can
confirm the centering from the intensities.
* 10.6 per-batch relative-B cross-validates on ASU-group parity, not
the frame parity the other surfaces use.
* 13 the resolution cutoff sits one reported-shell width PAST the
CC1/2 = 0.30 crossing, so data below 0.30 are kept.
* 7.1 the orientation-only prior penalises all three components of the
angle-axis vector.
* 10.1 rlp is the RECIPROCAL Lorentz factor, L = 1/rlp.
* 14.3 the model scaling is fitted over work and free reflections
alike, so R-free is free of refinement, not of the scaling fit.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Making azim_int_settings.high_q_recipA and
spot_finding_settings.high_resolution_limit optional changes the wire
format: when unset they are omitted rather than sent with a placeholder,
and a client generated from an older spec does j.at() on them. The
azimuthal-integration limit now defaults to unset, so a stock broker
omits it out of the box - the break needs no operator action to hit.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The first-pass spot cache constructed an MXAnalysisWithoutFPGA and an
AzimuthalIntegrationProfile inside the per-image lambda, so every cache
miss allocated a CUDA stream, the preprocessing buffer, the spot finder,
the azimuthal integrator and the Bragg engine, used them for one frame,
and freed them again - hundreds of times, serially, on the
--redo-rotation-spots path. Both worker loops already hoist the same
object out of their loop; only this path did not.
Build them once for the whole first pass, and only when spots actually
have to be found (with --reuse-rotation-spots there is nothing to
allocate). Reuse is safe because every azimuthal-integration path -
CPU, GPU and the fused adaptive engine - clears the caller's profile
before adding to it, so each frame's output is unchanged.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The expected-variance weights decompose an observation's sigma^2 into a
background part and a Poisson signal part, then rebuild the signal part
at the reflection's merged mean. The decomposition subtracted corr*I with
I taken as-is, so a negative I ADDED to the background part: an
observation at I = -1.5 with sigma^2 = 1 came out with a base variance of
2.7 rather than 1.
That inflates the variance of precisely the down-fluctuated observations
the correction exists for. Below about one photon they are then
under-weighted and the merged mean is biased high - the same direction of
error, in the same regime, that weighting by the observation's own sigma
produces. Subtract max(0, I) instead: a negative intensity has no Poisson
signal to remove.
Both users of the decomposition are fixed - the stills merge, where
expected-variance weighting is now the default, and the rotation combine
it was mirrored from, which had it first.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
