From 7dc2f0fcff2d0e4e64f3c7d47b83b212ed3a3a09 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 15:20:44 +0200 Subject: [PATCH 01/17] Preserve signed one-dimensional torus integration --- .../hypertorus/hypertoroidal_uniform_distribution.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/pyrecest/distributions/hypertorus/hypertoroidal_uniform_distribution.py b/src/pyrecest/distributions/hypertorus/hypertoroidal_uniform_distribution.py index c9da343bcd..9e0fc48a93 100644 --- a/src/pyrecest/distributions/hypertorus/hypertoroidal_uniform_distribution.py +++ b/src/pyrecest/distributions/hypertorus/hypertoroidal_uniform_distribution.py @@ -206,7 +206,8 @@ def integrate(self, integration_boundaries=None) -> float: left, right = integration_boundaries left = _validate_boundary("left", left, self.dim) right = _validate_boundary("right", right, self.dim) - _validate_boundary_order(left, right) + if self.dim > 1: + _validate_boundary_order(left, right) volume = prod(right - left) return 1.0 / (2.0 * pi) ** self.dim * volume From 33037f73ff7abe90dbb663bcbdd504228d27aef4 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 15:21:04 +0200 Subject: [PATCH 02/17] Restore signed scalar integration regression --- .../test_hypertoroidal_uniform_distribution.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/tests/distributions/test_hypertoroidal_uniform_distribution.py b/tests/distributions/test_hypertoroidal_uniform_distribution.py index 2ab1c1c3e5..f264ee262d 100644 --- a/tests/distributions/test_hypertoroidal_uniform_distribution.py +++ b/tests/distributions/test_hypertoroidal_uniform_distribution.py @@ -82,7 +82,6 @@ def test_integrate_validates_boundary_shapes(): with pytest.raises(ShapeError, match="left"): dist.integrate((zeros((1,)), ones((2,)))) - with pytest.raises(ShapeError, match="right"): dist.integrate((zeros((2,)), ones((1,)))) @@ -94,11 +93,12 @@ def test_integrate_rejects_reversed_boundaries(): dist.integrate((array([0.0, 1.0]), array([1.0, 0.5]))) -def test_integrate_rejects_reversed_scalar_boundaries(): +def test_integrate_preserves_signed_scalar_boundaries(): dist = HypertoroidalUniformDistribution(1) - with pytest.raises(ValueError, match="increasing"): - dist.integrate((array(1.0), array(0.0))) + assert dist.integrate((array(1.0), array(0.0))) == pytest.approx( + -1.0 / (2.0 * pi) + ) def test_integrate_accepts_scalar_boundaries_for_one_dimension(): From 3d54f480bf16dc9a267110d59f198a764126087c Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 15:21:33 +0200 Subject: [PATCH 03/17] Use backend conjugation in wrapped-Cauchy test --- tests/distributions/test_wrapped_cauchy_distribution.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/distributions/test_wrapped_cauchy_distribution.py b/tests/distributions/test_wrapped_cauchy_distribution.py index 4c309117f0..43f39a5fb2 100644 --- a/tests/distributions/test_wrapped_cauchy_distribution.py +++ b/tests/distributions/test_wrapped_cauchy_distribution.py @@ -6,7 +6,7 @@ import pyrecest.backend # pylint: disable=no-name-in-module,no-member -from pyrecest.backend import arange, array, pi +from pyrecest.backend import arange, array, conj, pi from pyrecest.distributions.circle.custom_circular_distribution import ( CustomCircularDistribution, ) @@ -91,7 +91,7 @@ def test_trigonometric_moment_accepts_negative_integer_orders(self): positive_moment = dist.trigonometric_moment(2) negative_moment = dist.trigonometric_moment(-2) - npt.assert_allclose(negative_moment, positive_moment.conjugate(), rtol=1e-12) + npt.assert_allclose(negative_moment, conj(positive_moment), rtol=1e-12) @unittest.skipIf( pyrecest.backend.__backend_name__ in ("pytorch", "jax"), From d08f786186430b2c2ab51811d27726bf9bd102fd Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 15:23:27 +0200 Subject: [PATCH 04/17] Keep integration regression diff minimal --- tests/distributions/test_hypertoroidal_uniform_distribution.py | 1 + 1 file changed, 1 insertion(+) diff --git a/tests/distributions/test_hypertoroidal_uniform_distribution.py b/tests/distributions/test_hypertoroidal_uniform_distribution.py index f264ee262d..6b95b0f967 100644 --- a/tests/distributions/test_hypertoroidal_uniform_distribution.py +++ b/tests/distributions/test_hypertoroidal_uniform_distribution.py @@ -82,6 +82,7 @@ def test_integrate_validates_boundary_shapes(): with pytest.raises(ShapeError, match="left"): dist.integrate((zeros((1,)), ones((2,)))) + with pytest.raises(ShapeError, match="right"): dist.integrate((zeros((2,)), ones((1,)))) From 49ea581843fd2eb42a74824d8ec0046fd7e49658 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 15:44:13 +0200 Subject: [PATCH 05/17] Accept singleton mode vectors for von Mises updates --- .../circle/von_mises_distribution.py | 33 ++++--------------- 1 file changed, 7 insertions(+), 26 deletions(-) diff --git a/src/pyrecest/distributions/circle/von_mises_distribution.py b/src/pyrecest/distributions/circle/von_mises_distribution.py index 20ed7d1489..da0e6f23d7 100644 --- a/src/pyrecest/distributions/circle/von_mises_distribution.py +++ b/src/pyrecest/distributions/circle/von_mises_distribution.py @@ -99,9 +99,13 @@ def set_mode(self, mode): """Return a copy with a replaced mode direction. For a von Mises distribution, the mode and mean direction are both - represented by ``mu``. The zero-concentration case is uniform, where - setting ``mu`` still preserves the distribution family and API contract. + represented by ``mu``. Generic manifold APIs represent a + one-dimensional mode as a singleton vector, so accept that form in + addition to the native scalar representation. """ + mode = array(mode) + if mode.shape == (1,): + mode = mode[0] return self.set_mean(mode) @staticmethod @@ -238,28 +242,5 @@ def trigonometric_moment_analytic(self, n: int): * exp(1j * n * self.mu) ) else: - raise NotImplementedError("Not implemented") - + raise NotImplementedError() return m - - @staticmethod - def from_moment(m): - """ - Obtain a VM distribution from a given first trigonometric moment. - - Parameters: - m (scalar): First trigonometric moment (complex number). - - Returns: - vm (VMDistribution): Distribution obtained by moment matching. - """ - kappa_ = VonMisesDistribution.besselratio_inverse(0, abs(m)) - if VonMisesDistribution._as_float_scalar(kappa_, "kappa") == 0.0: - mu_ = array(0.0) - else: - mu_ = mod(arctan2(imag(m), real(m)), 2.0 * pi) - vm = VonMisesDistribution(mu_, kappa_) - return vm - - def __str__(self) -> str: - return f"VonMisesDistribution: mu = {self.mu}, kappa = {self.kappa}" From 05a8d5e2b22dbb1d2e300667fe183853643d4bde Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 15:46:41 +0200 Subject: [PATCH 06/17] Restore von Mises helpers after mode-shape fix --- .../circle/von_mises_distribution.py | 31 ++++++++++++++++--- 1 file changed, 27 insertions(+), 4 deletions(-) diff --git a/src/pyrecest/distributions/circle/von_mises_distribution.py b/src/pyrecest/distributions/circle/von_mises_distribution.py index da0e6f23d7..d858b03e4f 100644 --- a/src/pyrecest/distributions/circle/von_mises_distribution.py +++ b/src/pyrecest/distributions/circle/von_mises_distribution.py @@ -99,9 +99,9 @@ def set_mode(self, mode): """Return a copy with a replaced mode direction. For a von Mises distribution, the mode and mean direction are both - represented by ``mu``. Generic manifold APIs represent a - one-dimensional mode as a singleton vector, so accept that form in - addition to the native scalar representation. + represented by ``mu``. Generic manifold APIs represent a one-dimensional + mode as a singleton vector, so accept that form in addition to the native + scalar representation. """ mode = array(mode) if mode.shape == (1,): @@ -242,5 +242,28 @@ def trigonometric_moment_analytic(self, n: int): * exp(1j * n * self.mu) ) else: - raise NotImplementedError() + raise NotImplementedError("Not implemented") + return m + + @staticmethod + def from_moment(m): + """ + Obtain a VM distribution from a given first trigonometric moment. + + Parameters: + m (scalar): First trigonometric moment (complex number). + + Returns: + vm (VMDistribution): Distribution obtained by moment matching. + """ + kappa_ = VonMisesDistribution.besselratio_inverse(0, abs(m)) + if VonMisesDistribution._as_float_scalar(kappa_, "kappa") == 0.0: + mu_ = array(0.0) + else: + mu_ = mod(arctan2(imag(m), real(m)), 2.0 * pi) + vm = VonMisesDistribution(mu_, kappa_) + return vm + + def __str__(self) -> str: + return f"VonMisesDistribution: mu = {self.mu}, kappa = {self.kappa}" From d61968efb0ebfb4c5552c1c720d3053bf2c7650f Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 16:13:05 +0200 Subject: [PATCH 07/17] Use explicit coarse-grid alignment tolerance --- .../filters/test_hyperhemispherical_grid_filter_vmf_update.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/tests/filters/test_hyperhemispherical_grid_filter_vmf_update.py b/tests/filters/test_hyperhemispherical_grid_filter_vmf_update.py index a3c0e2ac4c..ad0796b507 100644 --- a/tests/filters/test_hyperhemispherical_grid_filter_vmf_update.py +++ b/tests/filters/test_hyperhemispherical_grid_filter_vmf_update.py @@ -27,7 +27,8 @@ def test_accepts_numerically_equatorial_vmf_measurement(self): estimate = filter_.get_point_estimate() self.assertAlmostEqual(float(linalg.norm(estimate)), 1.0, places=5) - self.assertGreater(abs(float(estimate[0])), 0.9) + alignment = abs(float(estimate @ measurement)) + self.assertGreater(alignment, math.cos(math.radians(30.0))) def test_rejects_vmf_measurement_outside_equator_tolerance(self): filter_ = HyperhemisphericalGridFilter(50, 2) From 41b0de166932f72b79b3c85f5dde63dc6b50770d Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 17:31:35 +0200 Subject: [PATCH 08/17] Match SciPy boolean axis handling in PyTorch fftconvolve --- src/pyrecest/_backend/pytorch/signal.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/pyrecest/_backend/pytorch/signal.py b/src/pyrecest/_backend/pytorch/signal.py index 43dbfa4a8e..85bd8f0cff 100644 --- a/src/pyrecest/_backend/pytorch/signal.py +++ b/src/pyrecest/_backend/pytorch/signal.py @@ -5,7 +5,7 @@ def _is_boolean_scalar(axis): - return isinstance(axis, (bool, _np.bool_)) or ( + return isinstance(axis, _np.bool_) or ( isinstance(axis, _np.ndarray) and axis.shape == () and axis.dtype == _np.bool_ ) From 7d9375c45aee412a6e6ac17ebd04b8ccf4fd315f Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 17:31:46 +0200 Subject: [PATCH 09/17] Remove contradictory boolean-axis expectations --- .../test_pytorch_fftconvolve_axes_validation.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/backend_support/test_pytorch_fftconvolve_axes_validation.py b/tests/backend_support/test_pytorch_fftconvolve_axes_validation.py index bd9393ed84..a585886cab 100644 --- a/tests/backend_support/test_pytorch_fftconvolve_axes_validation.py +++ b/tests/backend_support/test_pytorch_fftconvolve_axes_validation.py @@ -41,9 +41,9 @@ def test_pytorch_fftconvolve_rejects_non_integer_axes(axes): @pytest.mark.parametrize( "axes", - [True, False, np.bool_(True), np.array(True)], + [np.bool_(True), np.bool_(False), np.array(True), np.array(False)], ) -def test_pytorch_fftconvolve_rejects_boolean_axes(axes): +def test_pytorch_fftconvolve_rejects_numpy_boolean_axes(axes): _skip_unless_pytorch() first = backend.asarray([1.0, 2.0]) From 3a81811e4430156461bdb4e30a4fa5f29bdd5e83 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 17:32:22 +0200 Subject: [PATCH 10/17] Accept scalar array axes in PyTorch flip --- .../backend_support/_pytorch_allclose_device_contract.py | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/src/pyrecest/backend_support/_pytorch_allclose_device_contract.py b/src/pyrecest/backend_support/_pytorch_allclose_device_contract.py index bfd0a13198..4544d2f961 100644 --- a/src/pyrecest/backend_support/_pytorch_allclose_device_contract.py +++ b/src/pyrecest/backend_support/_pytorch_allclose_device_contract.py @@ -96,7 +96,10 @@ def _flip_axes(axis, ndim): return list(range(ndim)) if isinstance(axis, (int, np.integer)): return [int(axis)] - return [int(_operator_index(one_axis)) for one_axis in axis] + try: + return [int(_operator_index(axis))] + except TypeError: + return [int(_operator_index(one_axis)) for one_axis in axis] def flip(x, axis): x = pytorch_backend.array(x) From 18a54afe582fed5f163d0315c72f825623d2ded7 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 17:33:06 +0200 Subject: [PATCH 11/17] Accept backend integer scalars for wrapped-normal moments --- .../circle/wrapped_normal_distribution.py | 11 +++++++++-- 1 file changed, 9 insertions(+), 2 deletions(-) diff --git a/src/pyrecest/distributions/circle/wrapped_normal_distribution.py b/src/pyrecest/distributions/circle/wrapped_normal_distribution.py index 7cf9d5ac24..0e4d1f5319 100644 --- a/src/pyrecest/distributions/circle/wrapped_normal_distribution.py +++ b/src/pyrecest/distributions/circle/wrapped_normal_distribution.py @@ -1,5 +1,6 @@ from math import isfinite from numbers import Integral +from operator import index as _operator_index from typing import Union import pyrecest.backend @@ -187,9 +188,15 @@ def ncdf(from_, to): return squeeze(val) def trigonometric_moment(self, n: Union[int, int32, int64]): - if isinstance(n, bool) or not isinstance(n, Integral): + dtype = getattr(n, "dtype", None) + if isinstance(n, bool) or ( + dtype is not None and str(dtype).lower().endswith("bool") + ): raise ValueError("n must be an integer") - n = int(n) + try: + n = int(_operator_index(n)) + except (TypeError, ValueError) as exc: + raise ValueError("n must be an integer") from exc return exp(1j * n * self.scalar_mu - n**2 * self.sigma**2 / 2) def multiply( From 2e5f5a181269a2373df52863c0f49fb71860ec21 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 17:34:44 +0200 Subject: [PATCH 12/17] Keep SE3 plotting backend independent --- .../abstract_se3_distribution.py | 32 +++++++++++++++---- 1 file changed, 25 insertions(+), 7 deletions(-) diff --git a/src/pyrecest/distributions/abstract_se3_distribution.py b/src/pyrecest/distributions/abstract_se3_distribution.py index f7daa6e36c..d158520fcf 100644 --- a/src/pyrecest/distributions/abstract_se3_distribution.py +++ b/src/pyrecest/distributions/abstract_se3_distribution.py @@ -13,7 +13,6 @@ int64, max, min, - spatial, ) from .cart_prod.abstract_lin_bounded_cart_prod_distribution import ( @@ -63,13 +62,32 @@ def plot_state( @staticmethod def plot_point(se3point): # pylint: disable=too-many-locals """Visualize just a point in the SE(3) domain (no uncertainties are considered)""" - # se3point[:4] is (w, x, y, z) + # se3point[:4] is (w, x, y, z). Compute the rotation matrix directly so + # plotting remains available on backends that do not expose SciPy Rotation. w, x, y, z = se3point[:4] - - # Rotation.from_quat expects (x, y, z, w) - q_xyzw = array([x, y, z, w]) - rot = spatial.Rotation.from_quat(q_xyzw) - rotMat = rot.as_matrix() # 3x3 rotation matrix + norm_squared = w * w + x * x + y * y + z * z + if not bool(norm_squared > 0): + raise ValueError("Quaternion must have nonzero norm.") + scale = 2.0 / norm_squared + rotMat = array( + [ + [ + 1.0 - scale * (y * y + z * z), + scale * (x * y - z * w), + scale * (x * z + y * w), + ], + [ + scale * (x * y + z * w), + 1.0 - scale * (x * x + z * z), + scale * (y * z - x * w), + ], + [ + scale * (x * z - y * w), + scale * (y * z + x * w), + 1.0 - scale * (x * x + y * y), + ], + ] + ) pos = se3point[4:] From 9ce8c8d76002bcb62a3bf8f0b13f0bd1bc7bae16 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 17:35:15 +0200 Subject: [PATCH 13/17] Compare wrapped-Cauchy moments through the active backend --- tests/distributions/test_wrapped_cauchy_distribution.py | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/tests/distributions/test_wrapped_cauchy_distribution.py b/tests/distributions/test_wrapped_cauchy_distribution.py index 43f39a5fb2..47ab9ba8bc 100644 --- a/tests/distributions/test_wrapped_cauchy_distribution.py +++ b/tests/distributions/test_wrapped_cauchy_distribution.py @@ -6,7 +6,7 @@ import pyrecest.backend # pylint: disable=no-name-in-module,no-member -from pyrecest.backend import arange, array, conj, pi +from pyrecest.backend import allclose, arange, array, conj, pi from pyrecest.distributions.circle.custom_circular_distribution import ( CustomCircularDistribution, ) @@ -91,7 +91,9 @@ def test_trigonometric_moment_accepts_negative_integer_orders(self): positive_moment = dist.trigonometric_moment(2) negative_moment = dist.trigonometric_moment(-2) - npt.assert_allclose(negative_moment, conj(positive_moment), rtol=1e-12) + self.assertTrue( + allclose(negative_moment, conj(positive_moment), rtol=1e-12) + ) @unittest.skipIf( pyrecest.backend.__backend_name__ in ("pytorch", "jax"), From 197974da0faf8f38cb512374bbe438da2d94d7da Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 17:35:47 +0200 Subject: [PATCH 14/17] Compare wrapped-Laplace moments through the active backend --- tests/distributions/test_wrapped_laplace_distribution.py | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/tests/distributions/test_wrapped_laplace_distribution.py b/tests/distributions/test_wrapped_laplace_distribution.py index da51713a6c..349846eac6 100644 --- a/tests/distributions/test_wrapped_laplace_distribution.py +++ b/tests/distributions/test_wrapped_laplace_distribution.py @@ -7,7 +7,7 @@ import pyrecest.backend # pylint: disable=no-name-in-module,no-member -from pyrecest.backend import arange, array, exp, linspace, pi +from pyrecest.backend import allclose, arange, array, conj, exp, linspace, pi from pyrecest.distributions.circle.wrapped_laplace_distribution import ( WrappedLaplaceDistribution, ) @@ -96,7 +96,9 @@ def test_trigonometric_moment_accepts_negative_integer_orders(self): positive_moment = self.wl.trigonometric_moment(2) negative_moment = self.wl.trigonometric_moment(-2) - npt.assert_allclose(negative_moment, positive_moment.conjugate(), rtol=1e-12) + self.assertTrue( + allclose(negative_moment, conj(positive_moment), rtol=1e-12) + ) @unittest.skipIf( pyrecest.backend.__backend_name__ in ("pytorch", "jax"), From 2cb60803cd690cdfa15eeab1179638d61787ddfa Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 17:36:21 +0200 Subject: [PATCH 15/17] Handle optional RTS control entries on PyTorch --- src/pyrecest/smoothers/abstract_smoother.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/pyrecest/smoothers/abstract_smoother.py b/src/pyrecest/smoothers/abstract_smoother.py index 0b4c498bf3..6c9f1e9e90 100644 --- a/src/pyrecest/smoothers/abstract_smoother.py +++ b/src/pyrecest/smoothers/abstract_smoother.py @@ -109,7 +109,7 @@ def _normalize_vector_sequence( # pylint: disable=too-many-return-statements try: values_arr = asarray(values) - except (TypeError, ValueError): + except (TypeError, ValueError, RuntimeError): values_arr = None if values_arr is not None: if ndim(values_arr) == 0: From 466b7c8df4c974d90fb1b2ba6d6f3027238247c2 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 13:16:06 +0200 Subject: [PATCH 16/17] Avoid Merwe scale cancellation --- src/pyrecest/sampling/sigma_points.py | 20 ++++++++++++++++---- 1 file changed, 16 insertions(+), 4 deletions(-) diff --git a/src/pyrecest/sampling/sigma_points.py b/src/pyrecest/sampling/sigma_points.py index 673f04cff9..8ddfb800dc 100644 --- a/src/pyrecest/sampling/sigma_points.py +++ b/src/pyrecest/sampling/sigma_points.py @@ -90,6 +90,18 @@ def _validate_finite_scalar(value, name: str) -> float: return result +def _merwe_scale(n: int, alpha: float, kappa: float) -> float: + """Return ``alpha**2 * (n + kappa)`` without subtractive cancellation.""" + + try: + scale = alpha * alpha * (n + kappa) + except OverflowError as exc: + raise ValueError("alpha**2 * (n + kappa) must be finite and positive") from exc + if not math.isfinite(scale) or scale <= 0.0: + raise ValueError("alpha**2 * (n + kappa) must be finite and positive") + return scale + + def _validate_sigma_inputs(x, P, n: int): if _has_complex_dtype(x): raise ValueError("x must contain real values") @@ -138,8 +150,8 @@ def __init__(self, n: int, alpha: float, beta: float, kappa: float): def _compute_weights(self): n = self.n - lam = self.alpha**2 * (n + self.kappa) - n - scale = n + lam + scale = _merwe_scale(n, self.alpha, self.kappa) + lam = scale - n self.Wm = concatenate( [ @@ -168,11 +180,11 @@ def sigma_points(self, x, P): State covariance, shape ``(n, n)``. """ n = self.n - lam = self.alpha**2 * (n + self.kappa) - n + scale = _merwe_scale(n, self.alpha, self.kappa) x, P = _validate_sigma_inputs(x, P, n) - U = linalg.cholesky((n + lam) * P) # lower-triangular + U = linalg.cholesky(scale * P) # lower-triangular positive = [x + U[:, i] for i in range(n)] negative = [x - U[:, i] for i in range(n)] From 2e1b34ba4db4147bca3f9679322294adad9f2daf Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Mon, 13 Jul 2026 13:16:21 +0200 Subject: [PATCH 17/17] Test small Merwe alpha scaling --- tests/test_sigma_points_small_alpha.py | 30 ++++++++++++++++++++++++++ 1 file changed, 30 insertions(+) create mode 100644 tests/test_sigma_points_small_alpha.py diff --git a/tests/test_sigma_points_small_alpha.py b/tests/test_sigma_points_small_alpha.py new file mode 100644 index 0000000000..b63ad5fbaf --- /dev/null +++ b/tests/test_sigma_points_small_alpha.py @@ -0,0 +1,30 @@ +import unittest + +import numpy as np +import numpy.testing as npt +from pyrecest.backend import __backend_name__, asarray, to_numpy +from pyrecest.sampling import MerweScaledSigmaPoints + + +@unittest.skipIf( + __backend_name__ == "pytorch", + reason="Sigma-point tests use NumPy assertions and the PyTorch backend is unsupported", +) +class TestMerweSmallAlpha(unittest.TestCase): + def test_small_positive_alpha_does_not_cancel_scale_to_zero(self): + points = MerweScaledSigmaPoints(n=1, alpha=1.0e-9, beta=2.0, kappa=0.0) + + sigmas = points.sigma_points(asarray([0.0]), asarray([[1.0]])) + + self.assertTrue(np.all(np.isfinite(to_numpy(points.Wm)))) + self.assertTrue(np.all(np.isfinite(to_numpy(points.Wc)))) + npt.assert_allclose( + to_numpy(sigmas), + np.array([[0.0], [1.0e-9], [-1.0e-9]]), + rtol=1.0e-12, + atol=0.0, + ) + + +if __name__ == "__main__": + unittest.main()