rugnux: fit the mosaicity from the strongest spots only

The per-image mosaicity MLE ran over the whole indexed spot list, so it rode
on --max-spots, which is an indexing budget. A spot is detected when
I_full * R(tau) clears the finder threshold, so selecting by intensity censors
on R(tau): a deeper list holds proportionally more large-|tau| partially
recorded spots and the fit widens with it. Raising the budget 250 -> 1000
widened sigma_M 0.059 -> 0.075 deg on a rotation dataset whose measured rocking
width says 0.054.

That is not cosmetic. An over-wide mosaicity mis-states every partiality in
scaling: forcing the mosaicity across that range moved the merge error model
from b 0.039 / ISa 26 to b 0.167 / ISa 6, and the space-group search lost a
genuine 422 with it, merging the crystal in 222 instead.

Cap the fit at the strongest 250 spots. FilterSpotsByCount leaves the list
strongest-first, so this selects exactly the spots a smaller --max-spots would,
and the mosaicity becomes invariant: 0.0538 deg at 250, 500, 1000 and 2000
spots, with the correct space group at each. Trimming or down-weighting the
tau tail does not work - the censoring is multiplicative in R(tau), so it
widens the whole distribution rather than adding a tail.

Battery over 37 crystals: exactly one change, the demoted crystal repaired
(33 space groups matching XDS -> 34). 23 of 37 are bit-identical, never
reaching 250 spots. Unaffected elsewhere: the default spot count is 250, and
stills have no goniometer so they return before the fit.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>

docs/CHANGELOG.md

@@ -9,6 +9,7 @@ This is an UNSTABLE release. It includes many experimental features, as well as
 * Bragg prediction: How far the predictor walks the lattice is a setting (`bragg_integration_settings.max_hkl`) instead of a fixed 100, derived per crystal offline from the refined cell (`--max-hkl` overrides); the broker keeps a fixed bootstrap so a live acquisition has a predictable per-image cost.
 * rugnux: De-novo **space-group search** substantially more robust - centering ranked by net absences and judged on absent-class strength, merohedral-twin over-promotion vetoed, and genuine high-symmetry groups recovered on weak data.
 * rugnux: The space-group search takes systematic absences from the merge of all observations, needs at least three control reflections on an axial row to claim a **screw axis**, and no longer alters the production merge.
+* rugnux: The per-image mosaicity is fitted from the strongest 250 spots, so the indexing spot budget no longer sets it.
 * rugnux: Stills **partiality post-refinement** added, on by default (`--simple-stills` disables); several non-helping stills scaling/detection knobs removed.
 * rugnux: **Scaling** hardened against a collapsed per-frame scale on the stills path as well as rotation, and the merged-sigma systematic floor is kept when ISa is too degenerate to report.
 * rugnux: `--min-image-cc` now works for rotation data (opt-in); new `--search-min-zeta` (rotation default 0.85); reports how close a symmetry axis lies to the spindle.

docs/CPU_DATA_ANALYSIS.md

@@ -743,6 +743,8 @@ For rotation data the mosaicity $\sigma_M$ is estimated by maximum likelihood fr
 
 The $\phi$ search window for the Bragg angle is set **wider than the oscillation**, so that reflections recorded at large rocking offset are included. These tail reflections carry most of the information about the mosaic width; a window limited to the oscillation range would truncate the $\tau$ distribution and bias $\sigma_M$ low.
 
+The fit uses only the **strongest 250 spots** of an image, whatever the indexing spot budget (`--max-spots`) is. A spot is detected when $I_\mathrm{full}R(\tau)$ clears the finder threshold, so selecting spots by intensity censors on $R(\tau)$: a deeper list holds proportionally more large-$\tau$ partially recorded spots and the fit widens with it. Left uncapped, $\sigma_M$ therefore tracks the spot budget rather than the crystal — and since an over-wide mosaicity mis-states every partiality, the merge degrades sharply with it.
+
 The estimated mosaicity feeds the rotation prediction (how many frames each reflection spans, §8.3) and the rotation partiality (§10.2). It is **held fixed during scaling**: in the per-image scale fit the mosaicity is degenerate with the scale $G$ (both rescale the predicted intensity), so refining it there is unstable. A correct mosaicity matters because it controls both how much of each rocking curve is captured and the partiality used to form fulls (§10.6); too small a value truncates the captured curve and over-peaks the partiality, degrading the combined fulls.
 
 ---

image_analysis/indexing/AnalyzeIndexing.cpp

@@ -237,12 +237,27 @@ namespace {
         const Coord S0 = experiment.GetScatteringVector();
         const float delta_phi_rad = deg_to_rad(axis.GetWedge_deg());
 
+        // Fit from the strongest spots only, never the whole list. A spot is detected when
+        // I_full * R(tau) clears the finder's threshold, so the deeper the spot list reaches the more
+        // large-|tau| partially-recorded spots it holds - and sigma_M is fitted from exactly that tau
+        // spread. The estimate therefore rides on the INDEXING budget (--max-spots): taking 1000 spots
+        // per image instead of 250 widened it 0.059 -> 0.075 deg on a rotation dataset whose measured
+        // rocking width says 0.054. Trimming or down-weighting the tau tail does not remove this - the
+        // selection is multiplicative in R(tau), so it widens the whole distribution, not just the tail.
+        // An over-wide mosaicity then mis-states every partiality in scaling, which is far from
+        // harmless: on that dataset the merge error model went b 0.039 -> 0.167 and ISa 26 -> 6, and
+        // the space-group search lost the true 422 with it. FilterSpotsByCount leaves the list
+        // strongest-first, so taking the head selects exactly the spots a smaller --max-spots would.
+        constexpr size_t MOSAICITY_FIT_SPOTS = 250;
+        const size_t n_fit = std::min(spots.size(), MOSAICITY_FIT_SPOTS);
+
         std::vector<float> tau_values;
         std::vector<float> zeta_values;
-        tau_values.reserve(spots.size());
-        zeta_values.reserve(spots.size());
+        tau_values.reserve(n_fit);
+        zeta_values.reserve(n_fit);
 
-        for (const auto &s : spots) {
+        for (size_t si = 0; si < n_fit; ++si) {
+            const auto &s = spots[si];
             if (!s.indexed)
                 continue;