A High Precision Metacalibration Response Correction

Overview

Metacalibration measures the response of an estimator to a small applied shear by literally applying it ([4, 12]). An observed image is deconvolved by the point spread function (PSF), sheared by $\pm$ a small step in each component, and reconvolved by a slightly larger round target PSF. Differencing a shear estimator $e$ from the sheared images gives the response matrix $\mathbf{R}$, and the calibrated shear is $g = \langle\mathbf{R} \rangle^{-1} \langle e \rangle$.

The difficulty addressed in this note is noise. The deconvolution divides by the PSF in Fourier space, which amplifies the high-$k$ noise; the shear then stretches that amplified noise anisotropically. The metacal’d image therefore has noise with a power $P(k)$ that is no longer isotropic but modulated as $\cos2\theta$ (spin-2), and the modulation differs between the metacal types ($\pm$ and noshear). Because the response is a spin-2 measurement, a type-dependent spin-2 noise anisotropy leaks into $\mathbf{R}$ and biases the recovered shear ([3, 12]).

Two ways to remove it are

fixnoise (exact, noisy). Take an independent noise field, rotate it by $90^\circ$, metacal it, rotate back and add it to the image ([12]). A $90^\circ$ rotation flips the sign of every spin-2 (and spin-6) anisotropy, so the sum is isotropic in those harmonics. But the added field doubles the noise power, $\sqrt{2}$ in the noise, $\times2$ in variance, in every mode, whether or not that mode needed filling.

fusion (this note). Add to each $k$-mode only the power needed to cancel the spin-2 effect. Most of the metacal noise anisotropy is small; fixnoise fully doubles the power at all $k$ to correct it. The fusion¹ correction measures how much power each metacal type lacks relative to the others (the deficit), and supplies just that, by filtering the rotated field down to the deficit instead of adding it whole. The response is corrected to the same level as fixnoise, but noise increase is $\sim1.02$ instead of $\sqrt{2}$ for a round PSF.²

The Flow

The structure of this note follows the algorithm flow, which is:

1. $k$-space metacal (§ 2): the operator $K_t$ for each METACAL type $t$ that deconvolves, shears, and reconvolves a field in GALSIM $k$-space [7].

2. The impulse transfer (§ 3): push an impulse through $K_t$ to produce the per-type noise power transfer $P_t = |K_t|^2$, and see the type-dependent $m=2/m=6$ modulation.

3. The harmonics that bias shear (§ 4): the response isolates the 90°-odd harmonics $m=2, 6$; $m=4$ is inert.

4. The common-level deficit (§ 5): the power $D_t$ that fills each type’s mean and $m=2/m=6$ modulation up to the common (max-over-types) level, so every type gets one isotropic level.

5. The rotated correction field (§ 6): the source of the added power in the sky frame.

6. The fusion filter (§7): $H_t = \sqrt{D_t/P_t^{\mathrm{rot}}}$, shaping the rotated field to the deficit.

7. Assembly and output noise (§§§ 7–§ 8): one shared $k$-grid, and the predicted variance increase $F_V$.

The end-to-end validation is in § 9.

Metacalibration in $k$-space

In this section we describe the basic process implemented in the metacal software package³. The same underlying operations are used in METADETECTION [11], the nickname for a METACALIBRATION that includes the object detection phase.

Write the observed image as $I = (G * P)$, the true surface brightness $G$ convolved by the PSF $P$. Define the metacal type $t \in\{\text{noshear}, 1p, 1m, 2p, 2m\}$ to correspond with a shear operator $S_t$ that is identity for noshear and a $\pm$ step shear in $g_1$ for 1p/1m and $g_2$ for 2p/2m. Then the METACAL’d image is

i.e. deconvolve the PSF, apply the shear in the world (sky) frame, and reconvolve by a common target $P_{\mathrm{tgt}}$. The target is supplied by the user⁴. Whatever the target, METACAL dilates it by $1 + 2$ step to suppress modes exposed by the shearing process.

All of (1) is done in GALSIM $k$-space. The image and PSF are wrapped as GALSIM InterpolatedImages (Lanczos-15 $x$-interpolant, the metacal kernel); the PSF is deconvolved, the result sheared, and the round target reconvolved; the single complex profile is drawn once with drawKImage onto a padded $N_p \times N_p$ grid for which the $k$-spacing, ifftshift reordering, and centering phase are chosen to match numpy’s FFT layout. This is followed by one inverse FFT and a crop to the central $N \times N$ to yield $I_t$. The pad size $N_p$ is GALSIM’s automatic drawFFT size; it is exposed only so that a set of METACAL’d images can be forced onto one shared grid (§ 7), which the correction requires

The noise transfer from an impulse

The METACALIBRATION operation (1) is linear in the input field, so in $k$-space it acts mode-by-mode as multiplication by a complex transfer $K_t(\mathbf{k})$:

If the input is a stationary noise field with white (flat) power spectrum, the output noise power is

the per-type noise power transfer.

The form of $K_t$ in (2) follows in two steps. First, a shear is a linear coordinate remap of the image plane, and under Fourier transform it stays a pointwise coordinate remap of $\mathbf{k}$:

the inverse shear acting on $\mathbf{k}$, with no amplitude factor because a pure shear is area-preserving with a unit Jacobian. In practice the shear is applied in the sky frame, so on the pixel $k$-grid the remap goes through the WCS distortion matrix, but it remains a unit-determinant pointwise remap. Because the remap acts pointwise it distributes over products and quotients, so Fourier-transforming (1) with $I * P^{-1} = \widehat{I}/\widehat{P}$ gives

The shear works on the numerator and the denominator. It must for a galaxy because shearing the galaxy is the point of METACAL. But a stationary white noise field has $\langle\lvert\widehat{n}(\mathbf{k})\rvert^2\rangle= \sigma^2$, flat in $\mathbf{k}$, so the remapped argument is unchanged. Sheared white noise is statistically identical to white noise at second order. Taking $\lvert\cdot\rvert^2$ of (5) therefore collapses all the shear dependence of the noise power onto the deconvolution kernel, giving (3) with $K_t$ as in (2). Shearing the deconvolved noise is statistically equivalent to deconvolving unsheared noise by the sheared PSF, $\widetilde{S_tP}(\mathbf{k}) = \widehat{P}(S_t^{-1}\mathbf{k})$. This is the sense in which the type-dependent anisotropy of $\mathcal{P}_t$ is created by the METACAL shear (Figure 1). It enters through the remapped argument $S_t^{-1}\mathbf{k}$, not through $\widehat{P}$, and so flips sign between 1p and 1m even for a round PSF.

The operator also includes the response of the Lanczos interpolant with which METACAL represents the image and the PSF (the parenthetical in (2)), and the padding, which is why we measure $K_t$ rather than evaluate the quotient analytically. This is distinct from any interpolation used to build a coadd: that acted on the noise before the METACAL operator, and is part of the input spectrum $N$, not $K_t$. See § 6.1 for the treatment of correlated noise.

We measure $\mathcal{P}_t$ directly by pushing a unit delta impulse through the same METACALIBRATION process as the data: a delta has $\lvert\widehat{\delta}(\mathbf{k})\rvert= 1$ in every mode, so the output $k$-array is $K_t$ and its squared modulus is $\mathcal{P}_t$. This is exact for the operator, including the Lanczos kernel and padding, not an analytic approximation. It is measured per METACAL type, capturing the small but crucial way the shear changes the noise power. The white-noise assumption is separate from the operator. For a correlated input with spectrum $N(\mathbf{k})$, equation (5) gives an output noise power $\lvert K_t(\mathbf{k})\rvert^2N(S_t^{-1}\mathbf{k})$, with the shear acting on $N$ as well (see § 6.1).

In Figure 1 we decompose the transfer for the 1p type (a $+g_1$ stretch) for a round PSF and a diagonal WCS. The transfer (top-left) is band-limited and looks round by eye, because the metacal step is only 0.01. Dividing out the isotropic envelope (top-right, the relative difference to noshear) exposes a spin-2 quadrupole. Power is enhanced along $k_x$ (the stretch axis) and suppressed along $k_y$, reaching $\pm10\%$ at the band edge. The formal $m=2$ and $m=6$ azimuthal projections (bottom) isolate the two harmonics the correction targets; $m=6$ is real but only $\sim5\%$ the size of $m=2$. This type-dependent $m=2$ exists even for a round PSF. It is created by the metacal shear, which is why it cannot be ignored.

Which harmonics bias the shear

Not every anisotropy of $\mathcal{P}_t$ matters. The metacal response in $g_1$ for estimator $e_1$ is

which is a finite difference of two METACAL shears that are $90^\circ$ rotations of each other (a $+g_1$ stretch along $k_x$ versus a $-g_1$ stretch along $k_y$). Hence $\mathbf{R}$ isolates the part of $\mathcal{P}_t$ that is odd under a $90^\circ$ rotation. A harmonic $\cos m\theta$ satisfies

i.e. $m=2,6,10,\ldots$ survive into $\mathbf{R}$, while the even-under-$90^\circ$ harmonics $m=0,4,8$ cancel. This is the reason the fusion fills $m=2$ and $m=6$ and leaves $m=4$ alone: those are the harmonics fixnoise removes (adding the true $90^\circ$-rotated field zeroes every $90^\circ$-odd harmonic), and $m=4$ is rotation-invariant and identical across types, so it does not bias the response directly.

The $m=6$ part is spin-6, formally orthogonal to a spin-2 estimator, yet it does bias the shear. The biasing quantity is the spin-2 moment of the weighted noise power, $\int\mathrm{d}^2 k\,W(\mathbf{k})e^{-2i\theta}\mathcal{P}(\mathbf{k})$, with $W$ the measurement kernel. Expanding both fields in harmonics, a noise harmonic $m$ produces a spin-2 only for a kernel harmonic $n=2-m$. An isotropic kernel ($n=0$) admits only $m=2$; but the square pixel grid and the Lanczos/fit weighting have an $m=4$ component ($n=\pm4$), which demodulates $m=6$ down to $m=2$ through the identity

So a spin-6 noise feature leaks into the shear by beating against a spin-4 anisotropy of the measurement. Filling $m=6$ removes that channel; an end-to-end test (§ 9) measures its effect at $8$–$9\sigma$.

The common-level harmonic deficit

We can only add noise power, so the correction will be a filtered noise field with sufficient power to cancel the problematic $m=2/6$ modes. The fields added to each METACAL type must have the same level in order for the response to be correct. Thus the goal is to raise every type up to a common level, even noshear which has no anisotropy for a round PSF. The common level is the elementwise maximum over types, annulus by annulus, of the azimuthal mean and of each filled $m$-mode amplitude. Filling each type up to that envelope makes the corrected transfers identical across types and results in isotropic total power in $m=2/m=6$.

For each radial annulus we project every type’s power onto its mean and its $(m=2,6)$ cosine/sine coefficients,

form the common envelopes over types,

and build the per-type deficit by (i) raising the mean to $\mu_\star$ and (ii) for each $m$, cancelling that type’s modulation and topping it up to the common amplitude $A_\star^{(m)}$ (here aligned to the type’s phase):

clipped at zero (a power deficit cannot be negative).

In Figure 2 we show the effect. The raw per-type transfers (top) have opposite quadrupoles. The 1p is stretched along $k_x$ and the 1m along $k_y$, with noshear nearly round. This is the type-dependent spin-2 that biases $\mathbf{R}$. Adding $D_t$ (bottom) puts all three types on one near-isotropic level; the faint residual ring is $m=4$ which does not produce a bias.

The azimuthal average $\langle\cdot\rangle$ needs to be calculated carefully in order to avoid a bias, which we discuss in Appendices § 13 and § 14.

The rotated correction field

The deficit $D_t$ gives how much power to add in each mode. The source of that power is a rotated, metacal’d noise field. Take an independent noise realization with the same statistics as the data, rotate it by $90^\circ$ in the sky frame, push it through the same metacal operations, and rotate the result back by $-90^\circ$. Its power transfer, measured by the same delta-impulse trick is the rotated transfer $\mathcal{P}^{\mathrm{rot}}_t(\mathbf{k}) = \lvert K_t(\mathcal{R}_{90}\mathbf{k})\rvert^2$ where $\mathcal{R}$ is the rotation operator.

A $90^\circ$ rotation maps $\cos2\theta\to-\cos2\theta$: it flips the sign of every spin-2 (and spin-6) modulation. So where the image noise transfer $\mathcal{P}_t$ has $m=2$ troughs (along $k_y$ for a $+g_1$ stretch), the rotated transfer $\mathcal{P}^{\mathrm{rot}}_t$ has $m=2$ peaks. These are the modes the deficit $D_t$ is used to fill. The rotated field is thus naturally aligned with the deficit. fixnoise adds all of it (doubling the variance), while the fusion keeps only the fraction that matches $D_t$. This alignment is visible in Figure 3 below, where $\mathcal{P}^{\mathrm{rot}}_{\mathrm{ns}}$ and $D_{\mathrm{ns}}$ share the same lobe orientation.

For a conformal WCS a sky $90^\circ$ and a pixel rotation (e.g. numpy.rot90) coincide. Under a non-conformal WCS they do not. A pixel-grid rotation is not a sky $90^\circ$, so it fails to flip the sky-frame spin-2 and leaves a spurious additive shear $c$. Rotating in the sky frame, before metacal, then rotating back removes it.

We now have, on one shared grid, the un-rotated transfer $\mathcal{P}_t$ (the noise to isotropize, with its deficit $D_t$) and the rotated transfer $\mathcal{P}^{\mathrm{rot}}_t$ (the correction field’s power, in the same final frame after the rotation back). The correction field is added through a filter $H_t(\mathbf{k})$ applied to its $k$-amplitude, so it contributes power $H_t^2\mathcal{P}^{\mathrm{rot}}_t$. Setting that equal to the deficit, $H_t^2\mathcal{P}^{\mathrm{rot}}_t=D_t$, gives

Adding $H_t\cdot$ (rotated field) then gives a total noise power $\mathcal{P}_t+H_t^2\mathcal{P}^{\mathrm{rot}}_t=\mathcal{P}_t+D_t$ on the common isotropic level of § 5.

In practice we need to be careful not to add more than the full rotated field in any mode, and roll $H$ smoothly to zero:

with $\epsilon=\mathrm{ktol}\cdot\max\mathcal{P}^{\mathrm{rot}}_t$; we find $\mathrm{ktol}=10^{-4}$ works well.

The cap $\min(\cdot,1)$; $H_t\leq1$, guarantees the fusion correction never adds more than the full rotated field in any mode. It can never add more noise than fixnoise. In practice the raw $\sqrt{D_t/\mathcal{P}^{\mathrm{rot}}_t}$ exceeds 1 only near the band edge (where $\mathcal{P}^{\mathrm{rot}}_t$ is collapsing), so for a round or mildly elliptical PSF the cap is essentially inactive in-band; it is a safety bound that engages only as the PSF becomes highly elliptical.

Out of band both $D_t$ and $\mathcal{P}^{\mathrm{rot}}_t$ go to zero, so $\sqrt{D_t/\mathcal{P}^{\mathrm{rot}}_t}$ is a noisy 0/0; the taper $\mathcal{P}^{\mathrm{rot}}_t/(\mathcal{P}^{\mathrm{rot}}_t+\epsilon)$ rolls $H_t$ smoothly to zero wherever $\mathcal{P}^{\mathrm{rot}}_t$ falls below $\epsilon$, preventing the filter from amplifying numerical noise where there is no signal.

In figure 3 we show the construction for an elliptical PSF with $g^{\mathrm{psf}}_1=0.01$ (reduced shear convention), comparable to the metacal step. The rotated transfer $\mathcal{P}^{\mathrm{rot}}_{\mathrm{ns}}$ (top-left) and the deficit $D_{\mathrm{ns}}$ (top-right) share the same lobe orientation, confirming the alignment argument of § 6; the filter $H_{\mathrm{ns}}$ (bottom-left) carves the rotated power down to the deficit. The $k$ cross-section (bottom-right) demonstrates the cap and taper. The raw $\sqrt{D/\mathcal{P}^{\mathrm{rot}}}$ rises toward the band edge but stays below the cap at 1. The cap engages only for a highly elliptical PSF (see Figure 4). The taper pulls the final $H$ to zero beyond the band.

In figure 4 we show the same construction for a highly elliptical PSF $g^{\mathrm{psf}}_1=0.15$. The deficit is now of order the rotated transfer, the raw $\sqrt{D/\mathcal{P}^{\mathrm{rot}}}$ reaches 1 over half the band, and the $\min(\cdot,1)$ saturates the filter at $H=1$. In those modes the fusion correction adds the full rotated field, approaching fixnoise. This is the regime behind the rising added noise shown in § 8.

In figure 5 we repeat the construction for the nominal case, a round PSF, shown for the 1p type, where the only anisotropy is the standard metacal step (0.01). The deficit is the two-lobed $m=2$ pattern along $k_y$, which is the trough the $+g_1$ stretch produces in $\mathcal{P}_{1p}$ (§ 3), topping out at $\sim7\%$ of the transfer peak. The filter is correspondingly near zero along the $k_x$ axis. No noise is added where the metacal shear already enhanced the power. On the $k_x=0$ cross-section through the lobes the $\sqrt{D/\mathcal{P}^{\mathrm{rot}}}$ stays below $\sim0.5$ everywhere, so the cap never engages; the added power $H^2\mathcal{P}^{\mathrm{rot}}$ is a few percent of the plain metacal noise. This is a geometrical explanation of the $\sim2\%$ round-PSF noise cost of § 8.

Correlated Input Noise

$\mathcal{P}_t$ and $\mathcal{P}^{\mathrm{rot}}_t$ are defined relative to white input noise, and $H_t$ is exact only for white noise. For a correlated input with spectrum $N(\mathbf{k})$ the image noise power is $\mathcal{P}_t(\mathbf{k})N(S_t^{-1}\mathbf{k})$ (§ 3). The $90^\circ$ rotation leaves an isotropic $N$ unchanged, but the METACAL shear does not. To first order in the shear $g_1$ of type $t$, the remap is $S_t^{-1}\mathbf{k}\approx(k_x(1-g_1),k_y(1+g_1))$, so for isotropic $N(k)$,

with $\phi$ the azimuth of $\mathbf{k}$. The remapped $N$ has a spin-2 modulation of fractional amplitude $g_1\,d\ln N/d\ln k$. Since $g_1=\pm$ step for 1p/1m (and the same with $\sin2\phi$ for 2p/2m), the modulation flips sign between the paired types. It is type-odd so enters the shear response (§ 4), and the white deficit $D_t$ does not include it. The second-order terms are type-common and cancel in $R$.

How much this matters depends on how much $N$ varies across the PSF band. For typical coadds, Lanczos-3 resampling at unit pixel ratio leaves $N$ flat to $\sim3\%$ across the band. A $k$-space calculation with the exact resampling kernel of a 10-epoch dithered coadd, with and without rotational dithers, gives a type-odd residual after the correction of $6\times10^{-4}$ of the noise power, compared to $4\times10^{-4}$ for white noise (1.2% vs. 0.8% of the uncorrected anisotropy), and based on focused simulations we expect a multiplicative bias of $m\sim10^{-3}$ for typical Lanczos-3 coadds (§ 9) below requirements for upcoming surveys (e.g. [5]).

For strongly correlated noise the white filter is not adequate. For white noise smoothed by a Gaussian with $\sigma=1$ pixel it leaves half the type-odd anisotropy in place (1.6$\times10^{-2}$ of $3.2 \times10^{-2}$). The remedy is to measure the transfers by pushing the noise correlation kernel $\sqrt{N}$ through the operator in place of the delta, after which the deficit machinery is exact again (residual $3 \times10^{-4}$ in the same calculation); the kernel must come from the known resampling, as the autocorrelation of a single noise realization is far too noisy. The real noise realization is still needed in any case as it accurately represents whatever non-stationary structure the data noise has (Appendix § 15).

The predicted $FV$ of § 8 weights the band by $N$; for correlated noise it is approximate.

Assembling the correction on one shared grid

The filter, the deficit, and the noise to be cancelled must all be in the same $k$-frame. The filter is applied in $k$ on the padded grid with no second FFT, so the de-aliased correction field is in the same frame as the constructed transfer.

So the final image for type $t$ is the plain metacal’d data plus the filtered, rotated, metacal’d noise:

each cropped to the central $N \times N$. $\mathcal{F}^{-1}$ represents the inverse Fourier transform⁵.

The output noise level

The correction adds anisotropy-filling power, so the corrected per-pixel noise is higher than plain metacal. By Parseval on the padded grid, the per-pixel variance of a filtered white field is the mean of its squared filter. The plain METACAL contributes $\langle P_t\rangle$ and the correction adds $\langle H_t^{2}\mathcal{P}_t^{\mathrm{rot}}\rangle$. The deficit gives every type a common noise level, and only the $m=0$ (azimuthal-mean) part sets the total variance (the $m=2,6$ fills integrate to zero over each annulus). Calculating from noshear,

Downstream this multiplies the noise variance used for uncertainty calculations⁶. It is 1 when no correction is applied.

In figure Figure 6 we plot the resulting noise (std) increase $\sqrt{F_V}$ versus PSF ellipticity in the reduced shear convention [2]. For a round PSF the variance factor is $\approx1.04$, which is only $\sim2\%$ in the noise, so the correction is nearly free. It rises toward fixnoise’s $\sqrt{2}$ in noise ($\times2$ in variance) only as the PSF becomes highly elliptical ($\sqrt{F_V} \approx1.34$ at $g_1^{\mathrm{psf}} = 0.2$). Deconvolving by an elliptical PSF makes the metacal noise more anisotropic, so more of the rotated field is needed to fill the deeper deficit, and the cost approaches the full doubling of fixnoise. The prediction (16) matches the measured ratio of corrected to plain noise variance to better than 5% in direct white-noise tests.

End-to-end validation

The pieces above were each validated in dedicated shear-bias simulations, based on those from Sheldon et al., [10] but with additional features to support simulated coadds with edges Sheldon, [9]. For the primary simulations, we used “edgy”-coadd noise with 10 input epochs Sheldon, [9], semi-realistic galaxy population [10, 8], a 1.1% non-conformal WCS plus a $10^\circ$ rotation, and a spatially variable, elliptical Moffat PSF Moffat, [6].

For the analysis we used and full metadetection with $k$-space metacal, sep detection Barbary et al., [1] and a maximum-likelihood Gaussian fit. We did not push a correlated noise kernel through the deficit machinery so, as noted in § 6.1, for coadd correlated noise from Lanczos-3 resampling we expect $m \sim0.001$. The results were

• Multiplicative bias $0.0004 < m < 0.0014$ at 99.7% confidence, within LSST tolerances ([5]) and consistent with the expectation for uncorrected correlated noise.

• In a paired, non-noise-cancelling run, we measured $\mathrm{var}^{\mathrm{fixnoise}}/\mathrm{var}^{\mathrm{fusion}} = 2.043 \pm0.005$: fixnoise costs essentially double the shear-noise variance of the fusion correction (the residual above 2 is the fusion’s $\sim2.6\%$ larger restored response, not extra image noise).

• A matched-seed A/B ($m=\{2,6\}$ vs. $m=2$-only, $1.25 \times10^{6}$ trials, 11% non-conformal WCS, sky-rotated correction) resolves $m=6$ at 8–9$\sigma$: dropping it leaves a $+3.5 \times10^{-4}$ multiplicative bias that the $m=6$ fill removes, consistent with the fixnoise reference.

Summary

The fusion noise correction makes the metacal noise spin-2 isotropic while adding the minimum amount of noise power. The flow is:

1. Push a unit delta impulse through the METACALIBRATION to get the per-type noise power transfer $\mathcal{P}_t = |K_t|^2$ (§ 3).

2. Form the common-level deficit $D_t$ that fills each METACAL type’s mean and $m=2/m=6$ modulation up to the max-over-types envelope.

3. Measure the rotated transfer $\mathcal{P}^{\mathrm{rot}}_t$ of a sky-90°-rotated, METACAL’d noise field (§ 6).

4. Filter that field by $H_t = \sqrt{D_t/\mathcal{P}^{\mathrm{rot}}_t}$ with a cap and a taper and add it to the METACAL’d data on one shared $k$-grid (§§7–§ 7).

The noise increase is $F_V = 1 + \langle H^2\mathcal{P}^{\mathrm{rot}}\rangle/\langle\mathcal{P}\rangle$, $\approx1.04$ for a round PSF and rising toward 2 only for highly elliptical PSFs (§ 8), compared to fixnoise’s $\times2$ ($\sqrt{2}$ in the noise). The corrected estimator is accurate to LSST requirements even under combined correlated non-stationary noise, non-conformal WCS, and PSF ellipticity.

Acknowledgments

We thank the computing staff of the Scientific Computing & Data Facility at Brookhaven National Laboratory for their support.

Appendices

The rotation-covariant azimuthal quadrature

On a square FFT grid an annulus is not a continuous circle but the finite set of lattice modes $|\boldsymbol{k}|$ that round into each radial bin. Those modes have the $D_4$ symmetry of the lattice. The lowest surviving anisotropic moment is therefore

with $Q > 0$ when modes cluster on the pixel axes, $Q < 0$ on the diagonals. This $Q$ bleeds into the discrete $m=2$ projection. Projecting a pure $m=2$ pattern of amplitude $A$ at orientation $\varphi$, $\mathcal{P}=A\cos2(\theta-\varphi)$, the integrand self-beats to a constant plus an $m=4$ term, and the $m=4$ term does not integrate to zero on the discrete ring because $\langle\cos4\theta\rangle=Q$:

The two METACAL shears make $m=2$ deficits with eigenaxes that are $45^\circ$ apart: $g_1$ (1p/1m) along the world axis, $g_2$ (2p/2m) along the diagonal. The biased projection therefore over-fills one and under-fills the other, leaving a response anisotropy $\delta=\frac{1}{2}(R_{11}-R_{22})\propto Q$. Rotating the WCS turns the world axes against the fixed grid, so $\delta\propto\cos4\theta_w$, causing a rotation leak.

We can replace the uniform discrete mean $\frac{1}{N}\sum_i$ by a weighted average so that the spurious moments are reduced. We chose an arc-length (trapezoidal / Voronoi) quadrature. Within each annulus, we sort the mode angles and weight each mode by the angular gap it occupies, $w_i\propto\frac{1}{2}(\theta_{i+1}-\theta_{i-1})$ cyclically, normalized to $\sum_i w_i=1$. This is the trapezoidal rule on the non-uniform angular grid. It approximates the continuum integral $\frac{1}{2\pi}\int\cos p\theta\,d\theta=0$ and results in $M_p=\langle\cos p\theta\rangle$ being reduced to the trapezoidal error for $m=4,8,\ldots$. It is pure grid geometry, so it is computed once per grid, with no noise.

In Figure 7 we demonstrate both the error and its removal. By injecting a unit $m=2$ at every orientation, the uniform mean ripples by $Q=0.090$, tracking $1+Q\cos4\varphi$ to the line, while the trapezoidal quadrature recovers a relatively flat amplitude with a $\sim25\times$ smaller ripple.

In validation simulations of § 9, we find a bias due to the $\cos4\theta_w$ leak of the uniform quadrature’s residual bias that changes sign with WCS angle, $-0.0012\rightarrow+0.0011$ ($\Delta\theta_w\approx+0.0024$, $\sim7\sigma$). The trapezoidal quadrature flattens it to $+0.0004\rightarrow+0.0004$ ($\Delta\theta_w\approx0$).

The sky-angle deficit projection

The second subtlety is which frame the harmonics are projected in. The shear bias is a sky-frame contraction $\operatorname{Tr}[qC]$ between the object’s quadrupole and the noise covariance, so it is the sky-frame spin-2 that must be cancelled, not the pixel-frame one. Under a non-conformal WCS (a pixel grid that is sheared relative to the sky) the two differ. For the projection we therefore bin and project in sky-frame angle. A pixel-frame mode $e^{i\mathbf{k}\cdot\mathbf{x}}$ with sky coordinates $\mathbf{u}=M\mathbf{x}$ (where $M$ is the pixel$\to$sky Jacobian) has sky wavevector $M^{-T}\mathbf{k}$; normalizing out the pixel scale (only the WCS shape matters),

and the annulus binning and the $\cos m\theta/\sin m\theta$ projections use $\mathbf{k}_{\mathrm{sky}}$. The Jacobian decomposes as $M=\mathrm{scale}\cdot R(\theta)\cdot S(g)$. A diagonal scale maps to the identity (a no-op) and a pure rotation merely rotates $\theta$, in which case the $(c_m,s_m)$ is a no-op. Only a genuine non-conformal shear changes the projection. Because the diagonal and pure-rotation cases reduce to identity automatically, the WCS Jacobian is all that is required. This sky-frame projection is the analogue of doing the metacal shear in the world frame, and of the sky-frame rotation of § 6. Every spin-2 operation in the code works in the sky frame.

The stationary-noise limit

The fusion filter uses only second moments of the noise, and the deficit $D_t$ is a non-negative power spectrum. If the data noise is stationary with a known spectrum $N(\mathbf{k})$, there is no need to filter a real realization. One can draw a single white field $\widetilde{w}(\mathbf{k})$ on the shared grid and add $\mathcal{F}^{-1}[\sqrt{D_t^N}\,\widetilde{w}]$ to each type, the same $\widetilde{w}$ for every type, where $D_t^N$ is the deficit of the correlated transfer $P_t(\mathbf{k})N(S_t^{-1}\mathbf{k})$, measured with the noise kernel $\sqrt{N}$ in place of the delta (§7); for white noise $D_t^N=D_t$. This removes the need for the noise image, the sky rotation, the second METACAL pass, $P_t^{\mathrm{rot}}$, the division and the cap. The variance increase is $1+\langle D_t^N\rangle/\langle P_tN\rangle$, which is $1+\langle D_t\rangle/\langle P_t\rangle$ for white noise. Where the cap is inactive the two constructions add the same power in every mode, so for Gaussian noise each type’s output noise is statistically identical; the correlation between types differs only through the phases of $K_t^{\mathrm{rot}}$, which are below 0.02 rad in band for the scenes below.

We compared them on the scenes of §7 with white Gaussian noise (2000 realizations, $48^2$ stamps). For the round PSF they agree. The predicted $m=2$ anisotropy of the 1p noise power falls from $4.8\times10^{-2}$ to $\approx5\times10^{-4}$ in both, the variance factor is 1.0428 in both, and the 1p−1m noise is the same. At $g_1^{\mathrm{psf}}=0.15$ the cap binds on 30–50% of the band and leaves 3–7% of the deficit unfilled, almost all of it $m=0$; both leave $\approx0.5\%$ of the original $m=2$. The Monte Carlo confirms both at the statistical floor ($\approx10^{-2}$). The synthesized field fills the common level completely, so its variance factor is 1.712 compared to the fusion’s 1.669, and in the stamp interior it adds correspondingly slightly more 1p−1m variance ($2.4\times10^{-3}$ compared to $\approx2\times10^{-3}$ per unit input variance). Over the full stamp the fusion’s 1p−1m variance is larger, because the sheared edge of the finite correction stamp enters at full strength where $H_t\approx1$; the synthesized field is stationary on the padded grid and has no edge.

This is a limiting case, not a replacement for the full fusion correction. Correlated noise is not necessarily the obstacle, since a known stationary $N$ enters through the kernel transfer. But $N$ must be known from the resampling that produced it; it cannot be measured well enough from a single realization. Furthermore, real coadd noise can be non-stationary (cell edges, varying depth, masks), in which case there is no single $N$ to synthesize from; the real realization represents that structure and $N$ (§7). The synthesized form is useful for simulations and as a check on the fusion machinery.