rotation: fix mosaicity underestimate + default-on capture-uncertainty (ISa 10.7->19.1)

The rotation per-image mosaicity was ~3x too small (0.045 vs the true 0.13deg),
which crippled the partiality model and capped per-observation precision: it
predicted reflections on too few frames and over-peaked the rocking partiality,
so the rot3d-combined fulls were ~1.7x noisier than XDS's, the integration
bottleneck behind the jungfraujoch-vs-XDS ISa gap.

Two root causes, both fixed:
- CalcMosaicityXDS (Kabsch-2010 MLE) searched each spot's exact-Bragg phi only
  within +-wedge (the 0.2deg oscillation). Reflections recorded at larger rocking
  offset - the tail that defines the mosaic width - fell outside and were dropped,
  truncating the tau distribution so the MLE underestimated ~2x. Widen the search
  window to wedge+0.8deg; the MLE then converges to the true 0.13deg (and is
  insensitive to widening further, since it weights by the recorded fraction).
- ScaleOnTheFly then re-refined the mosaicity from the intensity residual, which
  is degenerate with the per-image scale G and collapses it toward its floor.
  Trust the (now correct) indexing mosaicity and keep it fixed during scaling.

With the correct mosaicity, --capture-uncertainty (which down-weights the
over-extrapolated under-captured fulls) now pays off strongly, so default it ON
(1.0) for the rot3d combine; it stays off for non-rot3d. Together on the HEWL
rotation crystal: ISa 10.7 -> 19.1, and anomalous peak height vs XDS goes from
52% to ~78% (CL_CL 1.92x -> 1.29x). This reaches XDS's own published-correction
ceiling (DECAY+ABSORP+MODPIX ~= 19.6); the remaining gap to its quoted ISa 28 is
the I->inf extrapolation. No effect on the stills path (rotation-only code).

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>

image_analysis/indexing/AnalyzeIndexing.cpp

@@ -250,8 +250,14 @@ namespace {
                               + bstar * static_cast<float>(s.k)
                               + cstar * static_cast<float>(s.l);
 
-            // Find predicted phi angle
-            const auto phi_pred_opt = predict_phi_deg_local(pstar, S0, m2, 0.0f, axis.GetWedge_deg());
+            // Find predicted phi angle. The search window must be wide enough to catch reflections
+            // recorded at large rocking offset (|tau| up to ~mosaicity + dphi/2). Using ±wedge alone
+            // clips the tau tail at the oscillation width, so the MLE then underestimates the mosaicity
+            // ~2x (the tail reflections are exactly the ones that define the mosaic width). A generous
+            // window (oscillation + ~0.8deg rocking allowance) lets the tail in; the MLE is insensitive
+            // to making it wider still (it weights by the recorded fraction R(tau), which decays).
+            const float window_deg = axis.GetWedge_deg() + 0.8f;
+            const auto phi_pred_opt = predict_phi_deg_local(pstar, S0, m2, 0.0f, window_deg);
             if (!phi_pred_opt.has_value() || !std::isfinite(phi_pred_opt.value()))
                 continue;
 

image_analysis/scale_merge/ScaleOnTheFly.cpp

@@ -310,8 +310,10 @@ void ScaleOnTheFly::Scale(IntegrationOutcome &integration_outcome) const {
         } else {
             problem.SetParameterBlockConstant(&result.wedge);
         }
-        problem.SetParameterLowerBound(&result.mos, 0, s.GetMinMosaicity());
-        problem.SetParameterUpperBound(&result.mos, 0, s.GetMaxMosaicity());
+        // Trust the (now correctly estimated) indexing mosaicity; do NOT re-refine it here. The
+        // per-image scaling residual fit is degenerate between G and the mosaicity, and it collapses
+        // the mosaicity toward its floor (3x too small) which corrupts the partiality. Keep it fixed.
+        problem.SetParameterBlockConstant(&result.mos);
     }
 
     ceres::Solver::Options options;

tools/jfjoch_process.cpp

@@ -56,7 +56,7 @@ void print_usage() {
     std::cout << "   -w, --wedge[=num]               Refine image wedge during scaling with starting wedge value" << std::endl;
     std::cout << "   --scaling-high-resolution <num> High resolution limit for spot finding (default: no limit)" << std::endl;
     std::cout << "   --min-partiality <num>          Minimum partiality to accept reflection (default: 0.02)" << std::endl;
-    std::cout << "   --capture-uncertainty <num>     rot3d: systematic sigma ~num*(1-captured_fraction)*I on under-captured fulls (default: 0 = off)" << std::endl;
+    std::cout << "   --capture-uncertainty <num>     rot3d: systematic sigma ~num*(1-captured_fraction)*I on under-captured fulls (default: 1.0 for rot3d, 0 otherwise)" << std::endl;
     std::cout << "   --reject-outliers <num>         Per-observation merge outlier rejection, N sigma from the per-reflection median (default: off; e.g. 6, XDS/DIALS-style)" << std::endl;
     std::cout << "   --reject-delta-cchalf <num>     Per-crystal CC1/2-delta rejection: drop images with deltaCChalf below mean - N*stddev (default: off; e.g. 2.5)" << std::endl;
     std::cout << "   --min-image-cc <num>            Per-image CC limit in percent (default: no limit)" << std::endl;
@@ -288,7 +288,7 @@ int main(int argc, char **argv) {
     std::string ref_column;
     std::string dump_observations;   // diagnostic: dump unmerged -P rot3d fulls to this path
     double min_partiality = 0.02;
-    double capture_uncertainty_coeff = 0.0;
+    std::optional<double> capture_uncertainty_arg;   // explicit --capture-uncertainty; default depends on rot3d
     double min_image_cc = 0.0;
     int64_t scaling_iter = 3;
     std::optional<CrystalLattice> forced_rotation_lattice;
@@ -510,7 +510,7 @@ int main(int argc, char **argv) {
             min_partiality = std::stod(optarg);
             break;
         case OPT_CAPTURE_UNCERTAINTY:
-            capture_uncertainty_coeff = std::stod(optarg);
+            capture_uncertainty_arg = std::stod(optarg);
             break;
         case OPT_INTEGRATION_RADIUS:
             integration_radius_arg = optarg;
@@ -715,7 +715,10 @@ int main(int argc, char **argv) {
     if (wedge_for_scaling.has_value())
         scaling_settings.RotationWedgeForScaling(wedge_for_scaling);
     scaling_settings.MinPartiality(min_partiality);
-    scaling_settings.CaptureUncertaintyCoeff(capture_uncertainty_coeff);
+    // Capture-aware systematic sigma defaults ON (1.0) for the rot3d combine - it down-weights the
+    // over-extrapolated under-captured fulls and, with the mosaicity fix, lifts rotation ISa/anomalous
+    // substantially. Off for non-rot3d (no combine). An explicit --capture-uncertainty always wins.
+    scaling_settings.CaptureUncertaintyCoeff(capture_uncertainty_arg.value_or(combine_3d ? 1.0 : 0.0));
     scaling_settings.MinCCForImage(min_image_cc);
     if (outlier_reject_nsigma)
         scaling_settings.OutlierRejectNsigma(*outlier_reject_nsigma);