Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

gh-275: add tests for glass.lensing #454

Draft
wants to merge 8 commits into
base: main
Choose a base branch
from
Draft
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
101 changes: 72 additions & 29 deletions tests/test_lensing.py
Original file line number Diff line number Diff line change
Expand Up @@ -10,55 +10,55 @@
MultiPlaneConvergence,
RadialWindow,
deflect,
from_convergence,
multi_plane_matrix,
multi_plane_weights,
shear_from_convergence, # noqa: F401
)

if typing.TYPE_CHECKING:
from cosmology import Cosmology


@pytest.mark.parametrize("usecomplex", [True, False])
def test_deflect_nsew(usecomplex: bool) -> None: # noqa: FBT001
d = 5.0
r = np.radians(d)
def test_from_convergence(rng: np.random.Generator) -> None:
"""Add unit tests for :func:`from_convergence`."""
# l_max = 100 # noqa: ERA001
n_side = 32

def alpha(re: float, im: float, *, usecomplex: bool) -> complex | list[float]:
return re + 1j * im if usecomplex else [re, im]
# create a convergence map
kappa = rng.integers(10, size=healpix.nside2npix(n_side))

# north
lon, lat = deflect(0.0, 0.0, alpha(r, 0, usecomplex=usecomplex))
np.testing.assert_allclose([lon, lat], [0.0, d], atol=1e-15)
# check with all False

# south
lon, lat = deflect(0.0, 0.0, alpha(-r, 0, usecomplex=usecomplex))
np.testing.assert_allclose([lon, lat], [0.0, -d], atol=1e-15)
results = from_convergence(kappa)
np.testing.assert_array_equal(results, ())

# east
lon, lat = deflect(0.0, 0.0, alpha(0, r, usecomplex=usecomplex))
np.testing.assert_allclose([lon, lat], [-d, 0.0], atol=1e-15)
# check all combinations of potential, deflection, shear being True

# west
lon, lat = deflect(0.0, 0.0, alpha(0, -r, usecomplex=usecomplex))
np.testing.assert_allclose([lon, lat], [d, 0.0], atol=1e-15)
results = from_convergence(kappa, potential=True)
np.testing.assert_array_equal(len(results), 1)

results = from_convergence(kappa, deflection=True)
np.testing.assert_array_equal(len(results), 1)

def test_deflect_many(rng: np.random.Generator) -> None:
n = 1000
abs_alpha = rng.uniform(0, 2 * np.pi, size=n)
arg_alpha = rng.uniform(-np.pi, np.pi, size=n)
results = from_convergence(kappa, shear=True)
np.testing.assert_array_equal(len(results), 1)

lon_ = np.degrees(rng.uniform(-np.pi, np.pi, size=n))
lat_ = np.degrees(np.arcsin(rng.uniform(-1, 1, size=n)))
results = from_convergence(kappa, potential=True, deflection=True)
np.testing.assert_array_equal(len(results), 2)

lon, lat = deflect(lon_, lat_, abs_alpha * np.exp(1j * arg_alpha))
results = from_convergence(kappa, potential=True, shear=True)
np.testing.assert_array_equal(len(results), 2)

x_, y_, z_ = healpix.ang2vec(lon_, lat_, lonlat=True)
x, y, z = healpix.ang2vec(lon, lat, lonlat=True)
results = from_convergence(kappa, deflection=True, shear=True)
np.testing.assert_array_equal(len(results), 2)

dotp = x * x_ + y * y_ + z * z_
results = from_convergence(kappa, potential=True, deflection=True, shear=True)
np.testing.assert_array_equal(len(results), 3)

np.testing.assert_allclose(dotp, np.cos(abs_alpha))

def test_shear_from_convergence() -> None:
"""Add unit tests for :func:`shear_from_convergence`."""


def test_multi_plane_matrix(
Expand Down Expand Up @@ -107,3 +107,46 @@ def test_multi_plane_weights(
wmat = multi_plane_weights(weights, shells, cosmo)

np.testing.assert_allclose(np.einsum("ij,ik", wmat, deltas), kappa)


@pytest.mark.parametrize("usecomplex", [True, False])
def test_deflect_nsew(usecomplex: bool) -> None: # noqa: FBT001
d = 5.0
r = np.radians(d)

def alpha(re: float, im: float, *, usecomplex: bool) -> complex | list[float]:
return re + 1j * im if usecomplex else [re, im]

# north
lon, lat = deflect(0.0, 0.0, alpha(r, 0, usecomplex=usecomplex))
np.testing.assert_allclose([lon, lat], [0.0, d], atol=1e-15)

# south
lon, lat = deflect(0.0, 0.0, alpha(-r, 0, usecomplex=usecomplex))
np.testing.assert_allclose([lon, lat], [0.0, -d], atol=1e-15)

# east
lon, lat = deflect(0.0, 0.0, alpha(0, r, usecomplex=usecomplex))
np.testing.assert_allclose([lon, lat], [-d, 0.0], atol=1e-15)

# west
lon, lat = deflect(0.0, 0.0, alpha(0, -r, usecomplex=usecomplex))
np.testing.assert_allclose([lon, lat], [d, 0.0], atol=1e-15)


def test_deflect_many(rng: np.random.Generator) -> None:
n = 1000
abs_alpha = rng.uniform(0, 2 * np.pi, size=n)
arg_alpha = rng.uniform(-np.pi, np.pi, size=n)

lon_ = np.degrees(rng.uniform(-np.pi, np.pi, size=n))
lat_ = np.degrees(np.arcsin(rng.uniform(-1, 1, size=n)))

lon, lat = deflect(lon_, lat_, abs_alpha * np.exp(1j * arg_alpha))

x_, y_, z_ = healpix.ang2vec(lon_, lat_, lonlat=True)
x, y, z = healpix.ang2vec(lon, lat, lonlat=True)

dotp = x * x_ + y * y_ + z * z_

np.testing.assert_allclose(dotp, np.cos(abs_alpha))
3 changes: 2 additions & 1 deletion tests/test_observations.py
Original file line number Diff line number Diff line change
Expand Up @@ -9,6 +9,7 @@
tomo_nz_gausserr,
vmap_galactic_ecliptic,
)
import healpix


def test_vmap_galactic_ecliptic() -> None:
Expand All @@ -18,7 +19,7 @@ def test_vmap_galactic_ecliptic() -> None:
# check shape

vmap = vmap_galactic_ecliptic(n_side)
np.testing.assert_array_equal(len(vmap), 12 * n_side**2)
np.testing.assert_array_equal(len(vmap), healpix.nside2npix(n_side))

# no rotation

Expand Down
3 changes: 2 additions & 1 deletion tests/test_points.py
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,7 @@
positions_from_delta,
uniform_positions,
)
import healpix

if typing.TYPE_CHECKING:
import collections.abc
Expand Down Expand Up @@ -115,7 +116,7 @@ def test_loglinear_bias(rng: np.random.Generator) -> None:
def test_positions_from_delta(rng: np.random.Generator) -> None: # noqa: PLR0915
# create maps that saturate the batching in the function
nside = 128
npix = 12 * nside**2
npix = healpix.nside2npix(nside)

# case: single-dimensional input

Expand Down
Loading