Extra tests for healpy PR #1108: synalm Cl boundary check fix

test_pr1108_extra.py

import numpy as np
import healpy as hp
import sys

def test_synalm_cl_exactly_lmax_length():
    cl = np.ones(33, dtype=np.float64)
    lmax = 32
    np.random.seed(123)
    alm = hp.synalm(cl, lmax=lmax)
    assert alm.dtype == np.complex128

def test_synalm_cl_shorter_than_lmax():
    cl = np.ones(20, dtype=np.float64)
    lmax = 32
    np.random.seed(123)
    alm_short = hp.synalm(cl, lmax=lmax)
    cl_padded = np.concatenate([cl, np.zeros(lmax + 1 - len(cl))])
    np.random.seed(123)
    alm_padded = hp.synalm(cl_padded, lmax=lmax)
    np.testing.assert_allclose(alm_short, alm_padded, rtol=1e-12)

def test_synalm_cl_single_element():
    cl = np.array([1.0], dtype=np.float64)
    lmax = 16
    np.random.seed(456)
    alm = hp.synalm(cl, lmax=lmax)
    cl_padded = np.concatenate([cl, np.zeros(lmax)])
    np.random.seed(456)
    alm_padded = hp.synalm(cl_padded, lmax=lmax)
    np.testing.assert_allclose(alm, alm_padded, rtol=1e-12)

def test_synalm_multi_field_truncated():
    lmax = 16
    cl_tt = np.full(lmax + 1, 1e-5)
    cl_ee = np.full(lmax + 1, 2e-5)
    cl_bb = np.full(lmax + 1, 3e-5)
    cl_te = np.full(lmax + 1, 1e-6)
    cls_trunc = [cl_tt[:-1], cl_ee[:-1], cl_bb[:-1], cl_te[:-1]]
    cls_padded = [np.concatenate([c, [0.0]]) for c in cls_trunc]
    np.random.seed(101)
    alms_trunc = hp.synalm(cls_trunc, lmax=lmax)
    np.random.seed(101)
    alms_padded = hp.synalm(cls_padded, lmax=lmax)
    for i in range(len(alms_trunc)):
        np.testing.assert_allclose(alms_trunc[i], alms_padded[i], rtol=1e-12)

def test_synalm_different_truncation_per_field():
    lmax = 16
    cl_tt = np.full(lmax + 1, 1e-5)
    cl_ee = np.full(lmax + 1, 2e-5)
    cl_bb = np.full(lmax + 1, 3e-5)
    cl_te = np.full(lmax + 1, 1e-6)
    cls_trunc = [cl_tt[:10], cl_ee[:5], cl_bb[:15], cl_te[:2]]
    cls_padded = [
        np.concatenate([c, np.zeros(lmax + 1 - len(c))]) for c in cls_trunc
    ]
    np.random.seed(202)
    alms_trunc = hp.synalm(cls_trunc, lmax=lmax)
    np.random.seed(202)
    alms_padded = hp.synalm(cls_padded, lmax=lmax)
    for i in range(len(alms_trunc)):
        np.testing.assert_allclose(alms_trunc[i], alms_padded[i], rtol=1e-12)

def test_synfast_truncated_cl():
    nside = 32
    lmax = 64
    cl = np.full(lmax + 1, 1e-5, dtype=np.float64)
    cl_trunc = cl[:-1]
    np.random.seed(303)
    map_trunc = hp.synfast(cl_trunc, nside, lmax=lmax)
    cl_padded = np.concatenate([cl_trunc, [0.0]])
    np.random.seed(303)
    map_padded = hp.synfast(cl_padded, nside, lmax=lmax)
    np.testing.assert_allclose(map_trunc, map_padded, rtol=1e-12)

def test_synalm_high_lmax():
    lmax = 128
    cl = np.full(lmax, 1e-6, dtype=np.float64)
    np.random.seed(404)
    alm_trunc = hp.synalm(cl, lmax=lmax)
    cl_padded = np.concatenate([cl, [0.0]])
    np.random.seed(404)
    alm_padded = hp.synalm(cl_padded, lmax=lmax)
    np.testing.assert_allclose(alm_trunc, alm_padded, rtol=1e-12)

def test_synalm_float32_cl():
    lmax = 16
    cl = np.full(lmax + 1, 1e-5, dtype=np.float32)
    cl_trunc = cl[:-1]
    np.random.seed(505)
    alm_trunc = hp.synalm(cl_trunc, lmax=lmax)
    cl_padded = np.concatenate([cl_trunc, np.array([0.0], dtype=np.float32)])
    np.random.seed(505)
    alm_padded = hp.synalm(cl_padded, lmax=lmax)
    np.testing.assert_allclose(alm_trunc, alm_padded, rtol=1e-6)

def test_synalm_roundtrip():
    lmax = 32
    nside = 64
    cl_input = np.zeros(lmax + 1)
    cl_input[2:10] = 1e-3
    cl_trunc = cl_input[:-1]
    np.random.seed(707)
    alm = hp.synalm(cl_trunc, lmax=lmax)
    m = hp.alm2map(alm, nside)
    cl_out = hp.anafast(m, lmax=lmax)
    assert cl_out[0] < cl_out[2] * 0.01

def test_synalm_none_cl_entries():
    lmax = 16
    cl_trunc = np.array([1.0], dtype=np.float64)
    cls = [cl_trunc, None, None, None]
    np.random.seed(808)
    alms = hp.synalm(cls, lmax=lmax)
    assert len(alms) == 3

def test_synalm_sentinel_beyond_logical_end():
    lmax = 16
    cl_storage = np.full(lmax + 1, 1e-5, dtype=np.float64)
    cl_storage[lmax] = 1e6
    cl_truncated = cl_storage[:-1]
    np.random.seed(42)
    alm_truncated = hp.synalm(cl_truncated, lmax=lmax)
    np.random.seed(42)
    alm_zero_padded = hp.synalm(np.concatenate((cl_truncated, [0.0])), lmax=lmax)
    np.testing.assert_allclose(alm_truncated, alm_zero_padded)


if __name__ == "__main__":
    tests = [
        test_synalm_cl_exactly_lmax_length,
        test_synalm_cl_shorter_than_lmax,
        test_synalm_cl_single_element,
        test_synalm_multi_field_truncated,
        test_synalm_different_truncation_per_field,
        test_synfast_truncated_cl,
        test_synalm_high_lmax,
        test_synalm_float32_cl,
        test_synalm_roundtrip,
        test_synalm_none_cl_entries,
        test_synalm_sentinel_beyond_logical_end,
    ]
    passed = 0
    failed = 0
    for t in tests:
        try:
            t()
            print(f"PASS: {t.__name__}")
            passed += 1
        except Exception as e:
            print(f"FAIL: {t.__name__}: {e}")
            failed += 1
    print(f"\nResults: {passed} passed, {failed} failed out of {len(tests)} tests")
    sys.exit(1 if failed else 0)