An investigation into why hyperdrive peel leaves a large coherent residual on MWA EoR1 data,
and what actually fixes it. Short version: the field catalogue is fine, GLEAM-X DR2 confirms it,
hyperdrive's shapelet code is correct, and the Fornax A model over-predicts visibility amplitude
on the baselines MWA calibration uses by about 29%. Rescaling it cuts the coherent residual by a
quarter in amplitude, which is 44% in power.
| item | value |
|---|---|
| Field | MWA EoR1, phase centre RA 60°, Dec -30° |
| Band | centre channel 143, 167-197 MHz, 182.4 MHz centroid |
| Observations | 1443210064 (2025-09-29), 1452344160 (2026-01-13) |
| Pipeline | Birli → SSINS → hyperdrive di-calibrate → solutions-apply → hyperdrive peel |
| Sky model | srclist_pumav3_EoR0LoBES_EoR1pietro_CenA-GP_2023-11-07 |
| Calibration | -n 8000 --uvw-min 30l --max-iterations 300 --freq-average 40kHz |
The sky model uses LoBES (Lynch et al. 2021) for EoR0, Procopio et al. 2017 for the EoR1 field, and the Line et al. 2020 shapelet model for Fornax A.
On 1443210064, peel works well as an algorithm. Dirty-image RMS drops from 0.51 to 0.053 Jy/beam, and the 1000 ionospherically-subtracted sources subtract to about 1% of their flux.
But separating the residual into thermal and coherent parts, by frequency-differencing and by time-averaging, shows a floor that does not integrate down:
| quantity, |uv| > 50λ | value |
|---|---|
| calibrated sky, per visibility | 68.4 Jy |
| thermal noise, per visibility | 16.6 Jy |
| coherent residual floor | 6.2 Jy |
| floor as fraction of sky amplitude | 9% |
Peel's per-source gains have a median of 1.13, flat across three decades of source flux, across distance from the pointing centre, and across component type. Flatness across flux rules out confusion or CLEAN-bias. It is a flux-scale error. Peel absorbs it for the sources it ionospherically subtracts and leaves it for everything else.
The obvious suspect is the decade-old Procopio et al. 2017 catalogue used for EoR1. GLEAM-X DR2 (Ross et al. 2024, VizieR VIII/113) covers this field with 624,866 sources and in-band spectral fits, so it should be a straight upgrade.
Cross-matching 338 clean, isolated, single-component point sources brighter than 0.5 Jy:
| quantity | GLEAM-X vs puma | what the data wants |
|---|---|---|
| flux ratio | 1.039, scatter 0.041 | 1.135 |
| in-band (167-197 MHz) spectral index difference | +0.026 ± 0.003, scatter 0.049 | -0.31 |
Two independently-reduced catalogues agree to 4% in flux and 0.03 in spectral index, while the data disagrees with both by 13.5% and 0.31. The field catalogue is not the error. Swapping it closes about a third of the flux gap and none of the spectral gap.
Simulating Fornax A alone against the full 8000-source model gives its true share of
sum |model|², which is what weights the di-calibrate solution. A naive point-source estimate
from apparent flux says 92%, but Fornax A is heavily resolved, so the real number is much lower:
--uvw-min |
Fornax A weight | visibilities kept |
|---|---|---|
| 0λ | 61.0% | 100% |
| 30λ (operational) | 45.7% | 83.5% |
| 50λ | 35.6% | 72.8% |
| 100λ | 28.7% | 44.9% |
| 120λ | 19.6% | 26.0% |
| 150λ | 3.5% | 11.8% |
Still enough to set the flux scale. Least-squares fitting the model scale separately for Fornax A and for the field, after removing the other component from the data:
| observation | Fornax A sep / beam gain | Fornax A scale | field scale | required correction |
|---|---|---|---|---|
| 1443210064 | 13.5°, 0.248 | 0.924 | 1.160 | 0.796 |
| 1452344160 | 15.1°, 0.170 | 0.857 | 1.108 | 0.773 |
The required correction agrees to 3% across two nights while the beam gain differs by a factor of 1.46. That rules out a primary-beam error at large zenith angle. It is a source-model error, so one fix applies to every observation of this field.
Calibration splits the difference between a too-bright Fornax A and a correct field, which is why the field ends up ~13% under-subtracted and Fornax A slightly over-subtracted.
The error is scale-dependent in the same way on both nights, roughly 0.75 on 30-60λ rising to 0.90 beyond 100λ, so the model puts too much flux on the largest angular scales.
Fornax A is 7656 of the catalogue's 10280 shapelet coefficients, with basis orders up to n=86, while the next largest shapelet has 120 coefficients. A shapelet bug would therefore show up on Fornax A and nowhere else. Four checks, all clean:
- Basis table.
src/model/shapelet_basis_values.bin.gzis a 101 × 10001 f64 table. It matches the analytic Hermite-Gaussian to 1 part in 10¹⁴ for every order up to n=100. - Interpolation. Over the index range MWA actually touches, linear-interpolation error is 0.055% at n=86.
- Index clamping.
get_shapelet_envelopeinsrc/gpu/model.cuclamps negative indices to zero, a known hack added after Fornax A broke on SKA-Low-length baselines. For MWA the indices stay within 4858-5143 of a 0-10000 table, so it never triggers. - Normalisation and consistency. Extrapolating a no-beam Fornax A simulation to zero baseline
gives
|V| / S_catalogue = 0.9924. CPU and GPU agree to 7 parts in 10¹¹.
So the model genuinely over-predicts the visibility amplitude. The code reproduces it faithfully.
Six calibrations of 1452344160, each with --time-average 8s (which reproduces the operational 2s
solution to 2 parts in 10⁵ in mean gain while running in 50 minutes instead of 4 hours), then
solutions-apply, then compared against a vis-simulate model of the same source list. Coherent
residual is the non-thermal part on |uv| > 50λ.
| calibration | field flux-scale error | coherent residual | vs 30λ | gain scatter |
|---|---|---|---|---|
| 30λ, original model | 10.0% | 11.84 Jy | 0 | 1.067% |
| 50λ, original model | 7.2% | 11.28 Jy | -4.8% | 1.207% |
| 75λ, original model | 6.6% | 11.23 Jy | -5.2% | 1.272% |
| 120λ, original model | 5.3% | 11.49 Jy | -2.9% | 2.585%, 2 unconverged |
| 30λ, Fornax A ×0.773 | 0.7% | 8.90 Jy | -24.8% | 1.057% |
| 75λ, Fornax A ×0.773 | 0.0% | 8.80 Jy | -25.6% | 1.269% |
"Gain scatter" is the median channel-to-channel scatter of the solution amplitude, a proxy for solution noise.
Reading this:
- Raising
--uvw-mindoes work, and in the predicted direction: as Fornax A loses weight the field's flux-scale error falls from 10.0% to 5.3%, while Fornax A's own residual gets worse. - But the cutoff axis is only worth about 5%. It is flat from 50λ to 75λ and turns over by 120λ, where solution noise overtakes the benefit.
- Rescaling Fornax A is worth 25%, five times more, and it discards no baselines so solution smoothness is unchanged rather than degraded.
- Combining both adds only 0.8% over the rescale alone while costing 20% more gain scatter.
Scale the seven Fornax A components in the source list by about 0.77 and leave --uvw-min 30l
alone. This is a change to one source in the catalogue, costs nothing at runtime, and does not
change the component count, so the modelling budget is identical.
For a bulk reprocess this is the highest-value single change available: a 44% reduction in coherent residual power, for free.
- The correction is empirical and tuned to the baselines MWA samples. It deliberately makes the model's total flux disagree with Fornax A's integrated flux density. It fixes the visibilities you measure, not the source.
- Because the error is scale-dependent (0.75 short, 0.90 long), a flat rescale is a first-order fix. A proper refit of the Fornax A shapelet model would do better.
- After the fix the residual required correction is 0.754, so the optimum is nearer 0.752 than the 0.773 applied. One more iteration would gain a little.
- Verified on two observations of one field, consistent to 3%. Worth confirming on more before committing a bulk reprocess.
- Separately, there is a common-mode residual spectral tilt of about -0.3 in spectral index affecting the field and Fornax A alike. Per-source gains cannot absorb it because peel fits one gain per source for the whole 30 MHz band. This is unexplained and not addressed here.
# GLEAM-X DR2, whole catalogue as a 137 MB TSV in ~4 min.
# Column constraints are ignored, so you get all 624,866 rows.
curl -sS "https://cdsarc.cds.unistra.fr/viz-bin/asu-tsv?-source=VIII/113/catalog2\
&-out.max=unlimited&-out.add=_RAJ2000,_DEJ2000\
&-out=GLEAM-X,Fintwide,awide,bwide,psfawide,psfbwide,Fintfit200-SP,alpha-SP,chi2-SP,\
Fintfit200-CSP,alpha-CSP,beta-CSP,chi2-CSP" -o gleamx_dr2.tsv
# The CDS ftp catalog2.dat.gz is 891 MB and far slower.
# hyperdrive reads the puma FITS srclist identically to the YAML (200 MHz reference
# frequency), so edit the FITS: scale the FornaxA rows' NORM_COMP_PL / NORM_COMP_CPL.
# Fornax A's weight in calibration, measured rather than assumed:
hyperdrive vis-simulate -s srclist.fits --named-sources FornaxA -o fnxa.uvfits ...
hyperdrive vis-simulate -s srclist.fits -n 8000 -o all.uvfits ...
# then compare sum|V|^2 over the uv range di-calibrate actually uses.Gotcha that cost an hour: a wait loop using pgrep -f "hyperdrive di-calibrate" matches the
launching shell's own argv when the script is written via a heredoc, and deadlocks. Use
pgrep -x hyperdrive.