transform.resample_from_wcs() remakes gWCS if needed; supporting fixe…

…s required in astrodata.wcs

astrodata/wcs.py

@@ -104,11 +104,15 @@ def gwcs_to_fits(ndd, hdr=None):
     -------
     dict: values to insert into the FITS header to express this WCS
     """
+    if hdr is None:
+        hdr = {}
+
     wcs = ndd.wcs
     transform = wcs.forward_transform
     world_axes = list(wcs.output_frame.axes_names)
     nworld_axes = len(world_axes)
-    wcs_dict = {'WCSAXES': nworld_axes,
+    wcs_dict = {'NAXIS': len(ndd.shape),
+                'WCSAXES': nworld_axes,
                 'WCSDIM': nworld_axes}
     wcs_dict.update({f'CD{i+1}_{j+1}': 0. for j in range(nworld_axes)
                      for i in range(nworld_axes)})

gempy/library/transform.py

@@ -31,6 +31,7 @@
 from astropy.modeling.core import _model_oper, fix_inputs
 from astropy import table, units as u
 from astropy.wcs import WCS
+from astropy.io.fits import Header
 
 from gwcs import coordinate_frames as cf
 from gwcs.wcs import WCS as gWCS
@@ -43,6 +44,7 @@
 from geminidr.gemini.lookups import DQ_definitions as DQ
 
 import astrodata
+from astrodata import wcs as adwcs
 
 from .astromodels import Rotate2D, Shift2D, Scale2D
 from ..utils import logutils
@@ -527,7 +529,7 @@ def insert(self, index, model, copy=True):
             else:
                 self._models[i:i] = sequence
         # Update affinity based on new model (leave existing stuff alone)
-        self._affine &= astrodata.wcs.model_is_affine(model)
+        self._affine &= adwcs.model_is_affine(model)
 
     @staticmethod
     def split_compound_model(tree, ndim):
@@ -614,7 +616,7 @@ def __is_affine(self):
         """
         Test for affinity.
         """
-        return np.logical_and.reduce([astrodata.wcs.model_is_affine(m)
+        return np.logical_and.reduce([adwcs.model_is_affine(m)
                                       for m in self._models])
 
     @property
@@ -647,7 +649,7 @@ def affine_matrices(self, shape=None):
             except TypeError:  # self.ndim is None
                 raise TypeError("Cannot compute affine matrices without a "
                                 "dimensionality")
-        return astrodata.wcs.calculate_affine_matrices(self, shape)
+        return adwcs.calculate_affine_matrices(self, shape)
 
     def add_bounds(self, param, range):
         """Add bounds to a parameter"""
@@ -2060,6 +2062,16 @@ def resample_from_wcs(ad, frame_name, attributes=None, order=1, subsample=1,
             objcat['NUMBER'] = np.arange(len(objcat)) + 1
             setattr(ad_out[0], table_name, cat_table)
 
+    # We may need to remake the gWCS object. The issue here is with 2D spectra,
+    # where the resetting of the dispersion direction is done before the
+    # complete model and so isn't within the "WAVE" submodel. Converting to
+    # a FITS header and back results in a reconstructed model where "WAVE" is
+    # a distinct submodel.
+    wcs_fits_header = Header(adwcs.gwcs_to_fits(ad_out[0].nddata,
+                                                hdr=ad_out.phu))
+    if 'APPROXIMATE' not in wcs_fits_header.get('FITS-WCS', ''):
+        ad_out[0].wcs = adwcs.fitswcs_to_gwcs(wcs_fits_header)
+
     return ad_out