Skip to content

Covariance: credit Wu et al., fix chi_h^2, match cluster-lensing-cov, add validation ladder - #8

Merged
estevesjh merged 13 commits into
mainfrom
covariance-credits-and-chi-fix
Oct 8, 2026
Merged

estevesjh merged 13 commits into
mainfrom
covariance-credits-and-chi-fix

Conversation

@estevesjh

Copy link
Copy Markdown
Owner

Summary

Credits Hao-Yi Wu's paper and code, fixes a real bug in the Gaussian ΔΣ covariance, brings the shape noise in line with her latest commit, and adds a validation ladder against cluster-lensing-cov pinned at cddbb2a.

Credits (f2150a1): README, covariance docs and docs index now state that the covariance modules are a refactor and an optimized version of Wu et al. (2019, arXiv:1907.06611) and hywu/cluster-lensing-cov.

Bug fix (30f883f): DeltaSigmaGaussianCovariance integrated k dk instead of ℓ dℓ = χ_h² k dk, and the closed-form shot-shape term divided by the annulus area in Mpc² instead of the solid angle A/χ_h². Both were low by exactly 1/χ_h² (8e-7 at χ_h = 1100 Mpc) and invisible at χ_h = 1, where the old tests ran. Two tests that pinned the old behaviour were rewritten; a physical stack-variance test at χ_h = 1100 was added.

Shape noise (3ac9756): follows her cddbb2a convention — conditional ⟨Σcrit⟩/f_src with the 0.1 lens-source cut. Our old value was 12–40% too high (unnormalised ⟨Σcrit⟩ was divided by f_src again). f_src_behind gains min_separation/n_nodes; the kernels keep the 0.01 cut.

Optional Wu ℓ range (04c5ecd): ell_range="wu2019" reproduces her per-pair range [1/θ_max, 100/θ_min]. The converged integral stays the default, because her range drops about half of lss_lss/cross at the smallest r_p (none of it at ℓ < 10). A 12-bin covariance takes ~0.1 s (hers ~6.6 s).

Validation ladder (validation/): V0 snapshot of her code, V1 kernels, V2 Limber spectra, V3 noise, V4 Gaussian covariance term by term, V8 physical anchors, plus covariance_review_REPORT.md.

Results

  • 147 tests pass (test_covariance, test_fftlog_cov, test_limber, test_lensing_kernel, --no-cov).
  • V1 11 pass, V2 C_hh matches hers to ~1e-10, V3 pass (shape noise matches hers to 3e-16 once her rounded c = 3e5 km/s is removed), V4 20/20 pass (worst 1.2e-5), V8 pass.

Open decision (not changed here)

C_ΣΣ is ill-conditioned for foreground slabs: the source integral in q_Σ crosses a pole at z_s = z_h, so C_ΣΣ depends on slab width by up to 7.6e5×. Neither the paper nor her code says how sources in front of the halo should be treated. Two candidate definitions converge but differ by 3–7% on the largest-radius diagonal, where C_ΣΣ terms are 28–90% of the total. The default is unchanged; a question to Wu is drafted. Details in validation/covariance_review_REPORT.md.

Notes

  • The 400 kB reference snapshot (validation/data/clc_cddbb2a_reference.npz) is gitignored; regenerate with python validation/make_clc_reference.py (needs a cluster-lensing-cov clone, $CLUSTER_LENSING_COV_DIR; ~10 min, mostly her γ_t covariance). The ladder scripts need it.
  • Her repository has no license on GitHub; worth confirming with her that the derivation and credit wording are fine.
  • Not included: unrelated local edits to selection/scaling_relation.py and tests/test_selection.py.

🤖 Generated with Claude Code

estevesjh and others added 12 commits October 7, 2026 23:32
…dules

State that the covariance code is a refactor and an optimized version of
Hao-Yi Wu's method and code: acknowledgments and BibTeX in the README, a
provenance section in the covariance docs, and a credits section on the
docs index.

Co-Authored-By: Claude Sonnet 5.5 <noreply@anthropic.com>
The quadrature integrated k dk, but the covariance is an ell dell integral
and ell = k chi_h, so ell dell = chi_h^2 k dk. The noises are per
steradian, so the shot_shape closed form must divide by the solid angle
A_ann / chi_h^2, not the area in Mpc^2. Both were low by 1/chi_h^2
(8e-7 at chi_h = 1100 Mpc), invisible at chi_h = 1 where the tests ran.
This now matches cluster-lensing-cov, which integrates ell^2/2pi dln(ell).

Rewrite the two tests that pinned the old behaviour and add a physical
stack-variance test at a realistic distance, plus a quadrature-vs-closure
test at several chi_h. Correct the docstring and docs.

Co-Authored-By: Claude Sonnet 5.5 <noreply@anthropic.com>
…luster-lensing-cov cddbb2a)

LimberProjector.shape_noise_Sigma now forms
sigma_gamma^2/(n_s f_src) * (<Sigma_crit>/f_src)^2 with both integrals
starting at z_h + 0.1 (MIN_LENS_SOURCE_SEPARATION_NOISE), as her cddbb2a
does. The old form (unnormalised <Sigma_crit>, cut 0.01) was 12-40% high.
LensingKernel.f_src_behind gains min_separation/n_nodes (default unchanged);
the kernels q_sigma and <Sigma_crit^-1> keep the 0.01 cut.

V3 now agrees with her shape noise to 2e-14 once the exact c^4 and p(z)
normalisation constants are removed (raw -2.76e-3 = (c_ours/c_hers)^4 - 1).

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
…gma covariance

DeltaSigmaGaussianCovariance(ell_range="wu2019") reproduces her
_calc_C_ell_integration: per bin pair, ln ell on np.arange from
ln(1/theta_max) to ln(100/theta_min), step 1e-3, trapezoid. Pairs sharing
the outer bin share one grid, and each pair is a prefix of it, so the
integral is read off a cumulative sum (module-level
integrate_wu2019_pairs): ~0.1 s for 12 bins against her ~6.6 s.
The default stays the converged integral. Under "wu2019" shot_shape is
by default the quadrature over her range (her 0.9% truncation is part of
what is reproduced); exact_shot_shape=True overrides.

V4 now uses the public option: shot_shape agrees with hers to 1e-12, the
spectrum-carrying terms to 1.3e-5 (her stored spectra re-interpolated).

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
…in V2

For z_lss < z_h the signed q_Sigma source integral crosses z_s = z_h,
where Sigma_crit(z_s, z_h) has a simple pole: q_Sigma does not converge in
the node count and C_SigmaSigma moves by O(1) with DZ_SLAB. Restricting
sources to z_s > z_h removes both, but the paper and her code do not fix
which restriction (cut, normalisation) is intended, so defaults are
unchanged. Adds validation/diagnose_qsigma_pole.py (three candidate
definitions, their convergence and covariance effect), INFO rows in V2,
and short notes in the q_Sigma / limber docstrings.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
V1: her <Sigma_crit^-1> is exactly 0 at the last node (empty source
range), so the ratio was 0/0; compare nonzero nodes and require ours to be
exactly 0 elsewhere. V0: the growth scale-independence threshold becomes
1e-8 (the V2 C_hh tolerance it protects); measured 5.1e-9, below the
snapshot's 1.2e-7 interpolation error. Snapshot not regenerated.
README and docs/validation.md describe the covariance ladder and its
results.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Checks that do not use cluster-lensing-cov: shot x shape noise rebuilt from
the physical sample in arcmin units at several chi_h, a seeded Monte-Carlo
stack of shape noise, and the fraction of the line-of-sight integral that
lies below Wu et al.'s ell cut.

Co-Authored-By: Claude Sonnet 5.5 <noreply@anthropic.com>
f_src_behind, mean_sigma_crit and mean_inverse_sigma_crit integrated with a
100-node trapezoid. Against a 2048-node rule, 128 Gauss-Legendre nodes reach
~1e-12 on all three, while the trapezoid was off by 2e-4 (f_src), 5e-4 to
2e-3 (<Sigma_crit>, 0.1 cut) and 4-8% (<Sigma_crit>, 0.01 cut). The earlier
claim that <Sigma_crit> "does not converge when refined" was the trapezoid's
first interval, not the physics. q_sigma keeps the trapezoid: its integrand
has a pole at z_s = z_h for foreground slabs that no rule converges across.

Tests that pinned the trapezoid behaviour are rewritten (independent scipy
quad references, a convergence test, and a "converges but depends on the
cut" test replacing "does not converge").

Validation: V1, V3 and validate_lensing_kernel now check the formulas by
evaluating our integrands with her rule on her nodes (trapz_her_nodes, agrees
to floating point), and report the converged-vs-hers difference separately
as her quadrature error with a measured budget. The frozen-kernel tolerance
tightens from 3e-3 to 1e-6.

Co-Authored-By: Claude Sonnet 5.5 <noreply@anthropic.com>
Acknowledgments and BibTeX for McClintock et al. (2019), the DES Y3 cluster
cosmology paper (Abbott et al. 2025), cluster_toolkit and y3_cluster_cpp, in
the README and the docs credits.

Co-Authored-By: Claude Sonnet 5.5 <noreply@anthropic.com>
Use the references the y3_cluster_cpp docs give (the cluster count and
lensing forward model; the optical selection-bias and projection-lensing
model), point to github.com/estevesjh/y3_cluster_cpp, and drop the DES Y3
cosmology paper that was added without a source in those docs.

Co-Authored-By: Claude Sonnet 5.5 <noreply@anthropic.com>
96 nodes reach 3e-11 on <Sigma_crit> (0.01 cut, the worst case) and ~1e-12
or better on f_src and <Sigma_crit^-1> against a 2048-node rule; 128 was
more than needed. Updates the constant, its note, and the docs and report
that quoted the old number.

Co-Authored-By: Claude Sonnet 5.5 <noreply@anthropic.com>

@estevesjh estevesjh left a comment

Copy link
Copy Markdown
Owner Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Fixing trapz integrals; Adding missing citations

Comment thread src/clenspy/kernels/lensing_kernel.py Outdated
zs = self._zs_nodes(float(zh), min_separation, n_nodes)
if zs.size == 0:
continue
out[i] = np.trapezoid(self.survey.pz_src(zs), x=zs)

Copy link
Copy Markdown
Owner Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why we are using a trapz integration here? Please consider a more robust integration method, GL is a very good go to option for kernels.

Copy link
Copy Markdown
Owner Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Agreed, switched to Gauss–Legendre (8005552, node count set to 96 in 9fd3e55).

Measured against a 2048-node rule: 96 GL nodes reach 3e-11 on ⟨Σcrit⟩ (0.01 cut, the worst case) and ~1e-12 or better on f_src and ⟨Σcrit⁻¹⟩, while the 100-node trapezoid was off by 2e-4 (f_src), 5e-4–2e-3 (⟨Σcrit⟩, 0.1 cut) and 4–8% (⟨Σcrit⟩, 0.01 cut). f_src_behind, mean_sigma_crit and mean_inverse_sigma_crit now use gl_nodes (N_ZS_GL = 96).

q_sigma keeps the trapezoid: its integrand has a pole at z_s = z_h for foreground slabs that no rule converges across (the open decision in validation/covariance_review_REPORT.md).

Tests that pinned the trapezoid were rewritten against scipy quad, and V1/V3 now report the converged-vs-hers difference as her quadrature error.

Comment thread README.md
@@ -226,3 +161,33 @@ This project is licensed under the MIT License - see the [LICENSE](LICENSE) file
## Acknowledgments

CLensPy builds upon decades of research in weak gravitational lensing. We acknowledge the contributions of the broader weak lensing community to the theoretical foundations implemented in this package.

Copy link
Copy Markdown
Owner Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Add a citation to the DES cluster team, in particular for cluster tool-kit McClintock pacakge. And the y3-cluster-cpp github repo and paper.

Copy link
Copy Markdown
Owner Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Done in c632a7c and 38a731b. Acknowledgments and BibTeX (README) and the docs credits now cover McClintock et al. 2019 and cluster_toolkit, plus the two references the y3_cluster_cpp docs give: Aguena et al. 2023 (arXiv:2309.06593) and Costanzi et al. 2026 (PRD 113, 103508, arXiv:2604.05833), with the repo at github.com/estevesjh/y3_cluster_cpp.

Conflicts in docs/validation.md and validation/README.md were two sections
added in the same place by the two branches (the covariance ladder here, the
selection-bias / projected-density section on main); both are kept.

Co-Authored-By: Claude Sonnet 5.5 <noreply@anthropic.com>
@estevesjh
estevesjh merged commit 0afba18 into main Oct 8, 2026
9 checks passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant