Fix a bug preventing usage of numpy.linalg.inv and fix numpy 2.0 comp…

…atibility (#185)

* Fix a bug preventing usage of numpy.linalg.inv and fix numpy 2.0 compat

* Require a more recent version of sphinx_rtd_theme

* Update tweakwcs/linalg.py

Yes. Thanks!

Co-authored-by: Brett Graham <brettgraham@gmail.com>

* fix a bug in conversion to scalar

* fix regression test check for float epsilon

---------

Co-authored-by: Brett Graham <brettgraham@gmail.com>

CHANGELOG.rst

@@ -4,9 +4,19 @@
 Release Notes
 =============
 
-.. 0.8.3 (unreleased)
+.. 0.8.4 (unreleased)
    ==================
 
+0.8.3 (12-September-2023)
+=========================
+
+- Fixed a bug in the ``linalg`` module due to which computation of the inverse
+  matrix would fallback to custom implementation instead of using ``numpy``
+  implementation. [#185]
+
+- Fixed incompatibilities with the future (version 2.0) release of
+  ``numpy``. [#185]
+
 0.8.2 (13-April-2023)
 =====================
 

LICENSE.txt

@@ -1,4 +1,4 @@
-Copyright (C) 2018, Association of Universities for Research in Astronomy
+Copyright (C) 2023, Association of Universities for Research in Astronomy
 
 Redistribution and use in source and binary forms, with or without
 modification, are permitted provided that the following conditions are met:

docs/rtd-pip-requirements.txt

@@ -4,7 +4,7 @@ astropy
 numpydoc
 sphinx
 sphinx-automodapi
-sphinx-rtd-theme
+sphinx_rtd_theme>1.2.0
 stsci-rtd-theme
 sphinx-astropy
 graphviz

tweakwcs/correctors.py

@@ -474,9 +474,10 @@ def tanp_to_det(self, x, y):
             ra = np.atleast_1d(ra)
             dec = np.atleast_1d(dec)
         x, y = self._wcs.all_world2pix(ra, dec, 0, tolerance=1e-6, maxiter=50)
-        if not ndim:
-            x = float(x)
-            y = float(y)
+        if not ndim and np.size(x) == 1:
+            # convert to scalars:
+            x = float(x[0])
+            y = float(y[0])
         return x, y
 
     def world_to_tanp(self, ra, dec):

tweakwcs/linalg.py

@@ -56,7 +56,7 @@ def _is_longdouble_lte_flt_type(flt_type=np.double):
 
 def _find_max_linalg_type():
     max_type = None
-    for np_type in np.sctypes['float'][::-1]:  # pragma: no branch
+    for np_type in [np.longdouble, np.float64, np.float32]:  # pragma: no branch
         try:
             r = np.linalg.inv(np.identity(2, dtype=np_type))
             max_type = np_type
@@ -66,6 +66,7 @@ def _find_max_linalg_type():
             break
         except TypeError:
             continue
+    return max_type
 
 
 _USE_NUMPY_LINALG_INV = _is_longdouble_lte_flt_type(flt_type=np.double)
@@ -96,7 +97,7 @@ def inv(m):
     """
     # check that matrix is square:
     if _USE_NUMPY_LINALG_INV:
-        invm = np.linalg.inv(np.array(m, dtype=_MAX_LINALG_TYPE))
+        invm = np.linalg.inv(np.array(m).astype(_MAX_LINALG_TYPE))
         # detect singularity:
         if not np.all(np.isfinite(invm)):
             raise np.linalg.LinAlgError('Singular matrix.')

tweakwcs/tests/test_linearfit.py

@@ -6,6 +6,7 @@
 """
 from itertools import product
 import math
+import sys
 import pytest
 import numpy as np
 from astropy.modeling.models import Shift, Rotation2D
@@ -25,9 +26,12 @@
     'general': linearfit.fit_general,
 }
 
-_ATOL = 10 * _LARGE_SAMPLE_SIZE * np.sqrt(
-    np.finfo(linalg._MAX_LINALG_TYPE).eps
-)
+if linalg._MAX_LINALG_TYPE is not None:
+    _ATOL = 10 * _LARGE_SAMPLE_SIZE * np.sqrt(
+        np.finfo(linalg._MAX_LINALG_TYPE).eps
+    )
+else:
+    _ATOL = 10 * _LARGE_SAMPLE_SIZE * np.sqrt(sys.float_info.epsilon)
 
 
 def test_test_transform_selector():