Bugfix for the loss of array values

This was largely an issue of the zip(strict=False) causing missed entries. I made this strict=True

Author: WilliamJamieson

gwcs/coordinate_frames.py

@@ -388,6 +388,29 @@ def _native_world_axis_object_components(self):
         to be able to get the components in their native order.
         """
 
+    def add_units(self, arrays: u.Quantity | np.ndarray | float) -> tuple[u.Quantity]:
+        """
+        Add units to the arrays
+        """
+        return tuple(
+            u.Quantity(array, unit=unit)
+            for array, unit in zip(arrays, self.unit, strict=True)
+        )
+
+    def remove_units(
+        self, arrays: u.Quantity | np.ndarray | float
+    ) -> tuple[np.ndarray]:
+        """
+        Remove units from the input arrays
+        """
+        if self.naxes == 1:
+            arrays = (arrays,)
+
+        return tuple(
+            array.to_value(unit) if isinstance(array, u.Quantity) else array
+            for array, unit in zip(arrays, self.unit, strict=True)
+        )
+
 
 class CoordinateFrame(BaseCoordinateFrame):
     """

gwcs/wcs/_wcs.py

@@ -6,7 +6,6 @@
 
 import astropy.units as u
 import numpy as np
-from astropy import utils as astutil
 from astropy.io import fits
 from astropy.modeling import fix_inputs, projections
 from astropy.modeling.bounding_box import ModelBoundingBox as Bbox
@@ -117,24 +116,18 @@ def __init__(
         self._pixel_shape = None
 
     def _add_units_input(
-        self, arrays: list[np.ndarray], frame: CoordinateFrame | None
+        self, arrays: np.ndarray | float, frame: CoordinateFrame | None
     ) -> tuple[u.Quantity, ...]:
         if frame is not None:
-            return tuple(
-                u.Quantity(array, unit)
-                for array, unit in zip(arrays, frame.unit, strict=False)
-            )
+            return frame.add_units(arrays)
 
         return arrays
 
     def _remove_units_input(
         self, arrays: list[u.Quantity], frame: CoordinateFrame | None
     ) -> tuple[np.ndarray, ...]:
         if frame is not None:
-            return tuple(
-                array.to_value(unit) if isinstance(array, u.Quantity) else array
-                for array, unit in zip(arrays, frame.unit, strict=False)
-            )
+            return frame.remove_units(arrays)
 
         return arrays
 
@@ -166,10 +159,7 @@ def __call__(
         results = self._call_forward(
             *args, with_bounding_box=with_bounding_box, fill_value=fill_value, **kwargs
         )
-
         if with_units:
-            if not astutil.isiterable(results):
-                results = (results,)
             # values are always expected to be arrays or scalars not quantities
             results = self._remove_units_input(results, self.output_frame)
             high_level = values_to_high_level_objects(*results, low_level_wcs=self)
@@ -403,7 +393,9 @@ def outside_footprint(self, world_arrays):
                     max_ax = axis_range[~m].min()
                     outside = (coord > min_ax) & (coord < max_ax)
                 else:
-                    coord_ = self._remove_units_input([coord], self.output_frame)[0]
+                    coord_ = self._remove_quantity_output(
+                        world_arrays, self.output_frame
+                    )[idim]
                     outside = (coord_ < min_ax) | (coord_ > max_ax)
                 if np.any(outside):
                     if np.isscalar(coord):