Skip to content

Add Python test suite for qpoint package - #25

Open
arahlin wants to merge 7 commits into
fast-pix-exactfrom
add-test-suite
Open

arahlin wants to merge 7 commits into
fast-pix-exactfrom
add-test-suite

Conversation

@arahlin

@arahlin arahlin commented Sep 16, 2026 •

Copy link
Copy Markdown
Owner

qpoint had no tests. This adds a pytest suite — 448 of them — and the CI to run it, in six commits.

Based on fast-pix-exact (#34), and through it on fix-map-modes (#33).

No library behaviour changes here; src/qp_selftest.c is new C source but it is test
scaffolding. Faults the tests turned up were each fixed below as their own commit, on
master or on one of those two PRs — the tests for all of them are here, this being
where the suite lives:

  • Fix the map mode handling in init_dest and init_source #33, the map mode handling in init_dest/init_source, whose update=True path
    had no coverage at all and turned out to reinterpret the maps it was only meant to swap:
    TestUpdateDest, TestUpdateSource, TestVpolIsChecked and
    TestPolAndVpolDefaultToNone in test_qmap.py.
  • Make fast_pix exact, and turn it on by default #34, fast_pix, whose polarization angle lost 2e-10 near the poles and is now
    bit-identical to the angle path: TestFastPix in test_qpoint.py, which asserts that
    over the whole sky and inside each cap, and that the option is on by default. Eleven of
    its sixteen cases fail against the previous arithmetic, five of them only inside a cap.
tests
tests/test_qpoint.py 229 coordinate conversions, the parameter API, broadcasting, IERS Bulletin A, pixelization, galactic rotation
tests/test_qmap.py 185 mapmaking: flags, weights, gain, pair differencing, partial maps, derivative map modes, missing-pixel policy, threading, init_* map modes
tests/test_tools.py 8 the standalone refraction function
tests/test_astropy.py 11 azel2radec/radec2azel against astropy, sidereal time against its clock, the galactic rotation against ERFA
tests/test_selftest.py 2 one binding over src/qp_selftest.c
tests/test_benchmarks.py 13 timings, opt-in

Four of these pin behaviour that is easy to "fix" into something wrong:

  • Derivative map modes. Row count determines the mode and the mode determines which
    rows map2tod1 reads; the authority is the kernel, not qp_num_maps. The test perturbs
    every row of every valid row count and requires the TOD to change.
  • Threading. tod2map merges thread-local maps in arrival order, so only the hit
    counts come back identical. The rest is checked to a tolerance well inside the ~1e-15
    observed.
  • The solvers leave proj alone. None writes to a proj it does not hand back, which
    is what lets them read it in place; a write added on one of those paths would corrupt the
    caller's array silently.
  • IERS Bulletin A. The loader test uses three distinct constants, so a rotated column
    is visible rather than plausible.

Against astropy

test_astropy.py is the only external reference. The two agree to 14 mas once qpoint has
the same IERS data, which QPoint(update_iers=True) arranges for itself. The rates still
have to be on — rate_dut1/rate_wobble default to never, and the inverse transforms
need the _inv variants or radec2azel silently loses polar motion, a quarter of an
arcsecond. Without any of it the two disagree by 6 arcsec, which is also pinned. Nothing
touches the network: astropy's auto-download is switched off for the duration, and the
fixture is autouse because update_iers reads the table during construction.

Two transforms had no outside reference and now do. Sidereal time is checked against
astropy's, gmst and lmst both, longitude entering at lmst and nowhere else.

The galactic rotation is checked against erfa.icrs2g, which it matches exactly both
directions, rather than against astropy's Galactic, which it misses by 25 mas. That gap
is a convention difference — qpoint and ERFA use the Hipparcos angles straight from ICRS,
astropy reads the same three numbers as FK5 and takes the frame bias back on — and is
asserted separately so it stays on the record. A gal2radec round trip would catch none of
it, turning by the wrong angle and back still returning where it started.

The C constructors, exercised from C

qp_init_det, qp_init_detarr, qp_init_point, qp_init_map*, qp_copy_pixhash and
qp_copy_memory are public C API the Python package never calls — _libqpoint.py builds
those structs itself through its ctypes mirrors — so nothing a pytest can reach touches
them. src/qp_selftest.c drives them as an external C consumer would, and one binding runs
the lot: 59 numbered checks that each constructor fills the struct it claims to, that
copy=0 aliases and copy=1 copies, that a pixel hash survives being copied and still
refuses an absent pixel, that qp_copy_memory duplicates the Bulletin A table rather than
aliasing it, and that qp_reshape_map repoints every row when the row count changes. A
failure comes back as the number and text of the first check that failed.

Two constraints if you extend it: qp_init_map_pixhash wants the pixel list exactly as
long as the map, so the hash checks build a partial map; and the bulletin stores float,
so its test values are chosen to survive the narrowing.

Coverage

Measured by hand — coverage for the Python, gcov over a temporary -O0 --coverage
build for the C. All 400 pass under gcc with instrumentation.

qpoint C library
tools.py, _version.py 100% qp_selftest.c 100%
_libqpoint.py 93% qp_pixhash.c 97%
qpoint_class.py 90% qp_iers_bulletin_a.c 90%
qmap_class.py 87% qp_pixel.c 90%
total 90% qpoint.c 88%
qp_error.c 85%
qp_params.c 91%
qp_map.c 74%
quaternion.c 56%

The self-test is most of what moved in the C: qp_pixhash.c went 58% → 97% and
qp_iers_bulletin_a.c 79% → 90%, those being the files whose copy helpers are reachable
only through qp_copy_memory.

Nothing needed a gcov --exclude. qp_map.c read 47% before this suite, largely
because twenty-odd C-level constructors never ran. The two real cases were separated
instead: the ones the ctypes layer cannot reach are reached from C by the self-test, and
the five for which no declaration ever existed — qp_init_map_1d,
qp_init_map_from_arrays_1d and the three qp_init_detarr_*_from_array_1d — were deleted
on master, along with the four bulk setters, which is most of what moved qp_params.c.

quaternion.c's 56% is the remaining honest gap: declared, reachable public API with no
caller in the tree — the full 3×3 Quaternion_to_matrix, only _col3 being used, and the
QuaternionSlerp pair. An external C consumer can call those, so they count.

CI

.github/workflows/tests.yaml runs on ubuntu-latest and macos-latest against Python
3.9, 3.11 and 3.13. Optional dependencies are declared with markers and their tests skip
where absent, neither astropy nor healpy shipping wheels for 3.9:

test = ["pytest", "astropy>=6", "healpy>=1.17; python_version >= '3.10'"]

Every run prints the build it is testing, because it cannot be read off the result: the
threading tests pass on a serial build too, one thread trivially agreeing with itself, so a
build that quietly lost its parallelism would look exactly like one that kept it. There is
nothing in the C library to ask, so conftest.py reads the linkage of the built extension,
and on Linux under CI the suite fails if the build turns out serial. macOS is serial on
purpose: Apple clang ships no OpenMP runtime, and linking LLVM's collides with the one
healpy bundles.

Running it

pip install -e ".[test]"
pytest                  # ~1 s
pytest --benchmark      # adds the timing table

The benchmarks measure rather than assert, and a shared runner is too noisy for the numbers
to mean much. They are still collected without the flag, and each asserts the shape of its
result so one cannot quietly end up timing a no-op.

One detail in conftest.py that should not be undone: the timed call runs on a worker
thread. On the main thread, whose stack address is a deterministic function of the command
line, a tight loop's spilled temporaries can alias the heap buffers it streams — moving a
measurement by 3x, reproducibly enough to read as a result rather than as noise.

🤖 Generated with Claude Code

@arahlin arahlin self-assigned this Sep 16, 2026
@arahlin
arahlin requested a review from sbenton September 16, 2026 22:55
@arahlin
arahlin force-pushed the add-test-suite branch 14 times, most recently from 2e87cc7 to 31b6f29 Compare September 22, 2026 21:11
@arahlin
arahlin added this pull request to stack #28 September 23, 2026 16:53
@arahlin
arahlin force-pushed the add-test-suite branch 2 times, most recently from c099342 to 68ae844 Compare September 23, 2026 17:50
@arahlin
arahlin force-pushed the add-test-suite branch 9 times, most recently from 7f9da22 to 30db256 Compare September 24, 2026 02:33
arahlin and others added 7 commits September 30, 2026 12:18
qpoint had no tests. This is the suite over the shared API: coordinate
conversions and the quaternion round trips, the parameter API and its
grouping, IERS Bulletin A, pixelization, the galactic rotation, and
mapmaking through QMap -- flags, weights, gain, pair differencing,
partial maps, the derivative map modes and the missing-pixel policy.

Every per-sample argument is exercised as a scalar and as an array, in
both directions, because the two take different paths through the ctypes
layer: check_input broadcasts a scalar to full length while an array is
passed through, and the output shape then follows the longest argument.

.github/workflows/tests.yaml runs it on ubuntu-latest and macos-latest
against Python 3.9, 3.11 and 3.13. Optional dependencies are declared
with markers and their tests skip where absent, neither astropy nor
healpy publishing wheels for 3.9.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Five entry points and a set of mapmaking corners that the first pass
either skipped or exercised in one configuration only.

quat2pixpa carries its own copy of the fast_pix path, separate from the
one bore2pix takes, and nothing reached most of it: not the ordering it
picks itself, not the trig it picks, and not the branch where cos^2(b)
underflows near a pole and the angle has to come from the quaternion.
The exact poles go in to reach that branch, where only the angle is
required to be finite -- the pixel legitimately differs there, the angle
path's cos(theta) rounding to 1 and discarding the azimuth.

get_interp_val collapsed every map to the last one for a multi-map
input, and rotate_map raised on every call. load_bulletin_a could not
load a file at all: the columns came out rotated, and numpy's
unpack=True returns strided rows that ctypes would not take, so the
rotation was never reached. Its test uses three distinct constants, so a
rotated column is visible rather than plausible. Each fix is below this
branch, on master.

The mapmaking options are crossed rather than tried one at a time,
because the kernels branch on their combinations: a HWP alone and
differencing alone both miss the two bore2det_hwp calls inside
tod2map1_diff, a HWP on the binning side alone misses map2tod1's, and a
partial map reaches the pixel hash but not a failed lookup -- a skip
under error_missing=False and an error otherwise, decided again per
neighbour under interp_pix.

The depo had none of its own, though it is interface rather than
bookkeeping: every array QMap hands to C goes in there, which is also
what keeps the buffers alive, and spider_tools' UnifileMap and
MPIUnifileMap subclass QMap and read it. TestDepo pins the thirteen keys
and when each appears -- a switched-off component is False rather than
absent, the detector arrays exist only for the duration of a call, a
partial map carries its pixel list, each reset drops its own keys and
nothing else's, and an arbitrary key is allowed through so a subclass can
keep its own state there. TestDepoDownstreamPatterns replays what
spider_tools actually does with it, reduce_dest included, so the contract
it depends on is pinned here rather than only in its own repo.

get() with no arguments returns the four groups rather than a flat dict,
and each group is accepted back by set(). bore2radec refuses sindec
together with return_pa, the C having no entry point for it, and a call
with every argument scalar returns scalars.

Three behaviours are pinned here that nothing had stated. The first two
stand as they are; the third is a bug, fixed below this commit:

- from_tod tears the detector array down when it finishes, so the
  structure exists only between init_detarr and reset_detarr.
- polconv does not reach bore2radec. It negates beta in the kernel, so
  it shows up as the sign of the map's U row and nowhere else.
- solve_map_cho did not mask the pixels it cannot determine, returning
  -1509 for a two-hit pixel where the exact solver returned 0. With both
  solvers excluding the same pixels, the tests assert they agree
  everywhere and that both zero the rest.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Everything else in this suite checks qpoint against itself. These check
azel2radec and radec2azel against an independent implementation of the
same transforms, so a mistake has to be made twice to survive.

Matched properly the two agree to tens of milliarcseconds: 14 mas for
azel2radec, 43 mas for radec2azel, and 14 mas again with refraction
turned on for the same atmosphere. The tolerance is a tenth of an
arcsecond, which is loose enough not to chase rounding and tight enough
that a missing term cannot hide under it.

That 14 mas is a single term, not an accumulation, and the comment on
the tolerance now says which: astropy removes the Sun's gravitational
light deflection to recover the catalogue position and qpoint does not
model it at all, so qpoint's answer sits that far from the Sun. Over
800 pointings from 6 to 142 degrees elongation the residual projects
onto the Sun-radial direction with correlation 0.9997 and unit slope
against 4.07 mas * cot(elongation / 2), and its perpendicular component
averages 0.4 mas. Add the term back and the agreement is 0.8 mas, which
is the real precision of everything qpoint does model. Worth writing
down because the number otherwise reads as slack for rounding, and
because it grows toward the Sun: 4 mas at 90 degrees, 15 at 30, 78 at
6.

Matching them takes two steps that are worth writing down, since both
are easy to get wrong and quietly lose an arcsecond:

- qpoint needs the same Earth orientation data. The fixture writes a
  Bulletin A file out of the IERS-B table astropy ships with, so no test
  reaches the network, and rate_dut1 and rate_wobble have to be turned
  on for it to be used at all -- they default to never, and a loaded
  bulletin is otherwise read and ignored.
- The inverse transforms are switched separately, through rate_dut1_inv
  and rate_wobble_inv. Setting only the forward pair leaves radec2azel
  without polar motion, a quarter of an arcsecond, which a test asserts
  in both directions.

A last test measures what those terms are worth to a caller who never
loads a bulletin: arcseconds, not degrees.

Two more transforms had no outside reference of any kind, and get one
here. Sidereal time is checked against astropy's, gmst and lmst both:
longitude enters at lmst and nowhere else, so nothing else would catch
its sign. They need the Earth orientation data for the same reason the
coordinate transforms do -- sidereal time is Earth rotation, so
dropping dut1 costs the whole of that term, six arcseconds, which a
test pins as well.

The galactic rotation is checked against erfa.icrs2g rather than
astropy's Galactic frame, and agrees with it exactly, both ways. The
two references disagree by 25 mas and astropy is the one departing from
the canonical definition: the IAU 1958 system was defined against FK4
B1950, and Hipparcos supplies a rotation straight from ICRS with the
E-terms, the equinox and the ICRS-to-J2000 frame bias already folded
in, as three angles its catalogue calls exact. qpoint is built from
those angles and so is ERFA; astropy reads the same three numbers as
FK5 and goes ICRS -> FK5 -> Galactic, taking the frame bias back on.
The gap is asserted separately, so it stays on the record as a
convention difference and nobody corrects qpoint into matching astropy.
A round trip cannot see any of this: turning by the wrong angle and back
returns where it started.

astropy joins the test extra, so CI runs these rather than skipping
them. astropy 6 still supports every Python in the matrix, 3.9 included,
so every job gets the comparison.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
pytest --benchmark prints a comparison table; without the flag the
benchmarks are collected and skipped, which is enough to catch them
going stale. They measure rather than assert, a shared runner being too
noisy for the numbers to mean much, but each asserts the shape of its
result so one cannot quietly end up timing a no-op.

bore2pix is parametrised over fast_pix rather than left on the default,
so both paths are timed and the gap between the rows is what the option
is worth: 2.86 ms against 2.11 for 100k samples, about 1.35x.

The timed call runs on a worker thread, and has to stay there. The main
thread's stack sits below argv and the environment, so its address is a
deterministic function of the command line while every other mapping is
randomized, and at unlucky offsets a tight loop's spilled temporaries
alias the heap buffers it streams -- moving a measurement by 3x,
reproducibly enough to read as a result rather than as noise. A worker
thread gets a freshly mapped stack.

Every run prints the build it is testing, because that cannot be read
off the result: the threading tests pass on a serial build too, one
thread trivially agreeing with itself, so a build that quietly lost its
parallelism would look exactly like one that kept it. There is nothing
in the C library to ask, so conftest reads the linkage of the built
extension, and on Linux under CI the suite fails if the build turns out
serial. macOS is serial on purpose: Apple clang ships no OpenMP runtime,
and linking LLVM's collides with the one healpy bundles.

healpy joins the test extra for the pixelization comparisons, above
Python 3.9, where it publishes no wheels.

azel2bore and bore2azel are timed as a pair. They pay the same
corrections in opposite directions, and the inverse is where a rate
state that is not selected by direction shows up, so a row for one
without the other measures half of it. bore2azel closes the round trip
before reporting: its ctime comes last of four arguments where
azel2bore takes it last of seven, and passing it third is accepted in
silence -- the latitude arrives as ctime, so it is constant and
negative, every rate gate re-fires on every sample, and the row reads
19 us a sample instead of 106 ns. That looks like a library pathology
and is a caller error.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
qp_init_det, qp_init_detarr, qp_init_point, qp_init_map_from_arrays,
qp_copy_pixhash and qp_copy_memory are public C API the Python package
never calls: _libqpoint.py builds qp_det_t, qp_point_t and qp_map_t
itself, field by field through its ctypes mirrors. So nothing in the
suite reached them and their coverage read zero -- not because they were
known-good, but because they were unreachable from here.

src/qp_selftest.c drives them the way an external C consumer would, and
one binding runs the lot. Fifty-nine numbered checks: that the
constructors fill the structs they claim to, that copy=0 aliases the
caller's array and copy=1 copies it, that a pixel hash survives being
copied and still resolves and still rejects an absent pixel, that
qp_copy_memory duplicates the Bulletin A table rather than aliasing it
and the copy still interpolates, and that qp_reshape_map repoints every
row when num_vec or num_proj changes. Failures come back as the number
and text of the first check that failed, since a C abort would say
nothing.

qp_reshape_map is here for the same reason as the constructors: the 1-D
form it works on is reached only by assigning vec1d and its init flag
directly, which no constructor does and only the ctypes layer does, so
nothing a pytest can write touches it.

Check 28 covers qp_init_point marking ctime and q_hwp as malloc'd
whether or not it had allocated them, which left qp_free_point handing
free() two pointers the constructor never assigned -- undefined
behaviour on any call with time=0 or pol=0. The fix is below this
branch, on master.

Two constraints on extending it: qp_init_map_pixhash requires the pixel
list to be exactly as long as the map, so the hash checks build a partial
map rather than reusing the full-sky one; and the bulletin stores floats,
so its test values are chosen to survive the narrowing and keep the
comparison an equality.

qp_pixhash.c goes from 58% to 97% and qp_iers_bulletin_a.c from 79% to
90%, those being the files whose copy helpers are reachable only through
qp_copy_memory on a threaded build. Twenty-six functions in the family
now execute, sixteen of them fully.

Five more are defined in qp_map.c and declared nowhere -- qp_init_map_1d,
qp_init_map_from_arrays_1d and the three detarr *_from_array_1d variants.
No consumer can call them without inventing a declaration, so they are
left alone here: a deletion question rather than a testing one.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Coverage for five defects fixed on master and on the two branches under
this one, written here because tests/ only exists from here up. Each
group discriminates: it fails against the previous behaviour.

The inverse aberration rates did nothing. Both aberration corrections
read the forward rate state whichever direction they were called for, so
rate_aaber_inv and rate_daber_inv were inert, and because the check
writes ctime_last, a pass in one direction moved the other's clock. Two
tests fail against that; two more pass either way and stand as cover --
that an inverse pass cannot change what a forward pass then produces,
and that the round trip closes.

The rate-state sentinel was defeated by a negative ctime. 'once' has to
freeze whatever the epoch, but its "already done" guard was spelled
ctime_last > 0, which such a time never satisfies, so the state fell
through to the uninitialised branch and recomputed every sample --
indistinguishable from 'always'. The epoch is parametrised on purpose:
positive ctime passes either way, and so does a span straddling zero,
because once the time goes positive the guard starts firing. Three more
read the sentinel directly: fresh is NAN, a run leaves a real time
behind, reset_rates puts NAN back.

init_dest and init_source with update=True had no coverage at all, about
fifty lines including their shape checks. These take qmap_class.py from
77% to 85%, leaving only the reshape error branch, which needs a struct
no caller can build: the per-chunk zeroing and that accumulation
afterwards reproduces the same hit count, that the replacement reaches
the ctypes struct and not merely the depo, that init_dest returns
nothing, copy=True against the aliasing that is otherwise deliberate,
the partial-map path on both methods, a dest strict about its whole
shape and a source about npix alone, and that update on an uninitialized
QMap is an ordinary init.

One behaviour is pinned as it stands rather than as it should be: a
component switched off at init stays off and the map supplied for it is
dropped in silence.

pol and vpol default to None and are checked against the shapes rather
than used where the shapes settle the mode. Nine row/vpol combinations
on a source, the four contradictions init_dest refuses, that a supplied
proj stays authoritative where a supplied vec does not, and that pol
still decides at three rows.

fast_pix is exact, so TestFastPix asserts the polarization angle is
bit-identical to the angle path over the whole sky and inside each cap,
and the position angle likewise through quat2pixpa. Eleven of its
sixteen cases fail against the previous arithmetic, five only inside a
cap. It is on by default, for QPoint and QMap alike.

TestFastPixAgainstTheAnglePath takes the same claim through the public
entry points -- quat2pix, bore2pix, and quat2pixpa against the two steps
spelled out -- crossed with pol and the pixel ordering, since the fast
path picks the ordering itself. It also pins where the two are meant to
part: inside 2.1e-8 rad of a pole the angle path's cos(theta) rounds to
1 and dumps every direction into pixel 0, and the fast path is the one
that is right, so the test uses azimuths inside the quadrants rather
than on their boundaries, where the assignment is a tie-break that says
nothing. TestFastPixInMapmaking is the same check through the mapmaking
kernels, to_tod and from_tod.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
qp_update_ref reads sin(el) off the quaternion, which is unit only to
within rounding, so at el = 90 it could land a hair above 1 and asin
returned NaN. master clamps it now; this pins that for the forward
transform, where the NaN reached the boresight for 11664 of 172824
samples over a day of azimuth sweeps. The test fails with the clamp
removed.

The inverse never reproduced it -- it reads the elevation back after
undoing the aberration -- so it is left out rather than pinned by a
test that passes either way.

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

This branch has not been deployed

No deployments
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