reproject: geoid_height and itrf_transform silently accept mismatched lon/lat array shapes

[brendancol]

geoid_height and itrf_transform in xrspatial.reproject don't check that lon and lat arrays are the same shape. When they aren't, the numba parallel kernel reads past the end of the shorter array and the function returns wrong numbers without raising.

Reproducer

import numpy as np
from xrspatial.reproject import geoid_height, itrf_transform

# lon has 3 elements, lat has 2. Should raise; instead returns 3 values,
# the third using out-of-bounds lat memory.
print(geoid_height(np.array([0.0, 90.0, 45.0]),
                   np.array([0.0, 45.0])))
# [17.14536524 -59.66452408 -43.61396885]

print(itrf_transform(np.array([-74.0, 0.0, 45.0]),
                     np.array([40.7, 0.0]),
                     src='ITRF2014', tgt='ITRF2020', epoch=2024.0))
# 3-element result; the third entry is garbage

Why it matters

Both functions are public API. The numba kernels run under @njit(parallel=True), so OOB reads are undefined behaviour rather than an IndexError. The user gets plausible-looking numbers and no warning. You only catch it if you happen to write a test for it.

Proposed fix

Check lon.shape == lat.shape (and broadcastable h.shape) at the public-API boundary in geoid_height and itrf_transform, and raise ValueError with the mismatched shapes. The existing finite-value checks are the right place to add it.

Scope: _vertical.py (geoid_height) and _itrf.py (itrf_transform).

Tests: scalar+scalar still works, matching arrays don't regress, 1-D and 2-D mismatches raise, 0-D h still broadcasts against 1-D lon/lat.

Spotted during the reproject security sweep on 2026-05-17.