Stabilize point-set registration RMSE (#4965)

* Add one-shot registration RMSE fix

* Stabilize registration RMSE

---------

Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>

Author: FlorianPfaff

src/pyrecest/utils/_point_set_registration_common.py

@@ -9,7 +9,9 @@
 import numpy as np
 
 # pylint: disable=no-name-in-module,no-member
+from pyrecest.backend import abs as backend_abs
 from pyrecest.backend import all as backend_all
+from pyrecest.backend import max as backend_max
 from pyrecest.backend import (
     array_equal,
     asarray,
@@ -390,10 +392,15 @@ def default_cost(transformed_reference_points, moving_points):
 
 
 def compute_rmse(matched_costs) -> float:
-    """Compute the RMSE over matched costs."""
-    if matched_costs.shape[0] > 0:
-        return float(sqrt(mean(matched_costs * matched_costs)))
-    return float("inf")
+    """Compute the RMSE over matched costs without squaring overflow."""
+    if matched_costs.shape[0] == 0:
+        return float("inf")
+
+    scale = float(backend_max(backend_abs(matched_costs)))
+    if scale == 0.0 or not math.isfinite(scale):
+        return scale
+    scaled_costs = matched_costs / scale
+    return scale * float(sqrt(mean(scaled_costs * scaled_costs)))
 
 
 def summarize_assignment(assignment, costs) -> MatchSummary:

tests/test_point_set_registration_rmse_stability.py

@@ -0,0 +1,26 @@
+import math
+import unittest
+
+import numpy as np
+from pyrecest.backend import array, float64
+from pyrecest.utils._point_set_registration_common import compute_rmse
+
+
+class TestPointSetRegistrationRmseStability(unittest.TestCase):
+    def test_rmse_preserves_maximum_finite_cost(self):
+        largest = np.finfo(np.float64).max
+
+        rmse = compute_rmse(array([largest], dtype=float64))
+
+        self.assertTrue(math.isfinite(rmse))
+        self.assertEqual(rmse, largest)
+
+    def test_rmse_preserves_zero_and_empty_contract(self):
+        self.assertEqual(compute_rmse(array([0.0], dtype=float64)), 0.0)
+        self.assertTrue(
+            math.isinf(compute_rmse(array([], dtype=float64)))
+        )
+
+
+if __name__ == "__main__":
+    unittest.main()