Weight Q/U by their variance after the Stokes I division, and blank pixels with no Stokes I model

[AlecThomson]

Summary

Dividing Q/U by a Stokes I model divides their noise too. In channel j the corrected data have noise sigma_j * I(nu_0) / I_j, but the transform kept 1/sigma_j^2 weights. Channels where the model is faint got magnified noise at full weight, so PI blew up wherever a model fell towards zero. Only the reported noise (fractional_theoretical_noise) knew about the division.

stokes_i_weighting makes the noise-based weights (variance, natural, uniform_lsq, briggs) follow the division:

• "global" (default): w_j = T_j^2 / sigma_j^2 with T = nu^alpha, the same for every pixel. The map keeps one RMSF and one effective frequency. stokes_i_weight_alpha pins alpha, or "auto" fits it to the field's mean Stokes I (one extra pass over the Stokes I cube).
• "per_pixel": w_j = I_j^2 / sigma_j^2 with each pixel's own model. Inverse-variance weighting of the corrected data, so the most SNR, but the RMSF and effective frequency follow the spectral index. Forces the per-pixel RMSF, and cannot be combined with lam_sq_0_m2="per_pixel".
• None: the old weights.

Pixels with no usable Stokes I model are blanked (3D). Below stokes_i_snr_cut, with no Stokes I error, with Stokes I blank, or with a rejected model, the pixel used to divide by a flat model, i.e. no correction. It then reported a band-averaged PI next to fitted pixels reporting PI at the reference frequency, and the maps stepped at the edge of the fitted region. On the SB56289 cutout, with the cut raised to put that edge on the emission, the step was about 10%. Those pixels now get a NaN model, so their FDF and every map made from it are NaN, as in RM-Tools (rmtools_fitIcube leaves unfitted pixels NaN). Run without stokes_i for PI at every pixel; flint writes that set as fdf_no_i. fallback_model and StokesIFitOptions.fallback_alpha are gone.

Other changes:

• Blank pixels are NaN in the dirty and clean maps. calc_faraday_moments returns NaN for an all-NaN FDF rather than mom0 = 0 (a spectrum the threshold empties still gives 0). The 3D RM-CLEAN model stays at zero there, since a blank pixel has no components. compute_theoretical_noise gives NaN for a pixel with no weight, without a divide warning. These also apply to NaN mosaic edges.
• lambda^2_0. lam_sq_0_m2="auto" is the weighted mean of lambda^2 under the template weights in both modes, so REFFREQ moves.
• Noise. It still comes from the true Q/U sigma. noise_weights takes the noise weights separately from the transform weights.
• Results. RMSynth3DResults gains stokes_i_weighting and stokes_i_weight_alpha.
• 1D tool. run_rmsynth takes the same weighting options. A single spectrum has no map to step across, so 1D keeps its flat fallback for an unusable model.
• Default. The new default changes outputs for any run with a Stokes I model and a noise-based weight type.

Why

On SB56289 the very high PI came from order-2 Stokes I fits falling to a few percent of the reference flux at 1.8 GHz. On a 49x49 real-data cutout at fit_order=2, pixels with model gain above 10 went from 2.99x the order-1 PI to 1.52x (global) and 1.10x (per_pixel). Max PI fell from 4.48 to 2.38 and 2.25 mJy, and median SNR rose from 6.5 to 8.0.

Tests

• tests/test_stokes_i_weighting.py, new, covers:
  • option validation;
  • the theoretical noise against a Monte-Carlo in each mode;
  • SNR ordering None < global < per_pixel on a spectrum that fades;
  • lambda^2_0 from the template weights;
  • one RMSF and one FDF shape in global mode, and alpha-dependent RMSFs in per-pixel mode;
  • pixels below the SNR cut are NaN in every dirty and clean map and have a zero CLEAN model, in each mode, and fitted pixels report the right fraction;
  • field alpha recovery;
  • 1D/3D agreement under variance weighting in both modes;
  • one Stokes I fit per chunk in per-pixel mode with RM-CLEAN.
• The fallback tests in tests/test_tools_3d_stokes_i.py and tests/test_tools_3d_dask.py now check blanking instead: below the cut, no Stokes I error, blank or negative Stokes I, a runaway fit, a rejected model. Each checks that neighbouring pixels keep their fit.
• On the cutout, fitted pixels are bit-for-bit unchanged by the blanking, in dirty, clean and model peak PI.
• uv run pytest on the touched modules: 550 passed. The Stokes I notebook passes. uv run prek run --hook-stage manual --all-files passes.

Docs

docs/examples/rmsynth_3d_stokes_i.ipynb:

• "Masking low-SNR pixels" and "Blanked Stokes I" now say those pixels are blanked, and why, and show the uncorrected run covering them.
• "Models that run away" shows the rejected pixel blanked.
• New section "Weights after the Stokes I division", on a mock like the SB56289 band: the input, the SNR of the recovered PI and the RMSF width in each mode, and a table of the tradeoffs.

Please review the prose in those sections.

🤖 Generated with Claude Code

https://claude.ai/code/session_01F6ajgfHAviQKM4UXZ3gNi7

[codecov]

Codecov Report

❌ Patch coverage is 95.93496% with 5 lines in your changes missing coverage. Please review.
✅ Project coverage is 92.51%. Comparing base (d846352) to head (dfabc25).
✅ All tests successful. No failed tests found.

Files with missing lines	Patch %	Lines
rm_lite/tools_1d/rmsynth.py	92.85%	2 Missing ⚠️
rm_lite/utils/fitting.py	93.33%	2 Missing ⚠️
rm_lite/tools_3d/rmsynth.py	97.22%	1 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main      #97      +/-   ##
==========================================
+ Coverage   92.38%   92.51%   +0.12%     
==========================================
  Files          14       14              
  Lines        3114     3207      +93     
==========================================
+ Hits         2877     2967      +90     
- Misses        237      240       +3     

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

[AlecThomson]

Check Python 3.12 on ubuntu-latest failed on dfabc25 in test_rmclean_graph_build_stays_linear_in_chunk_count, a wall-clock test. It asserts that building the RM-CLEAN graph over 4x the chunks takes under 5x as long. The runner measured 5.46x (0.063 s vs 0.346 s).

It is not this PR's: rm_lite/tools_3d/rmclean.py is identical at dfabc25 and 724feb2, where this test passed in CI, and the maps change only runs inside the one map_blocks task per block, so it adds nothing to the graph. Locally, over 10 trials each, the ratio is 3.5 to 4.5 at dfabc25 and 3.5 to 4.9 at 724feb2. There is no fix to port. I have re-run the failed job once.

---

Generated by Claude Code