Merge branch 'master' into enh/pattern_noise_efficiency

# Conflicts:
#	geminidr/gnirs/recipes/sq/recipes_LS_SPECT.py

.jenkins/scripts/setup_agent.sh

@@ -18,4 +18,4 @@ fi
 source .jenkins/scripts/download_and_install_anaconda.sh
 
 conda install --yes pip wheel
-pip install "tox<=3.17.1" tox-conda
+pip install "tox>=3.8.1" tox-conda

astrodata/core.py

@@ -6,7 +6,7 @@
 import warnings
 from collections import OrderedDict
 from contextlib import suppress
-from copy import deepcopy
+from copy import copy, deepcopy
 from functools import partial
 
 import numpy as np
@@ -98,6 +98,25 @@ def __init__(self, nddata=None, tables=None, phu=None, indices=None,
         self._orig_filename = None
         self._path = None
 
+    def __copy__(self):
+        """
+        Returns a shallow copy of this instance.
+        """
+        obj = self.__class__()
+
+        for attr in ('_phu', '_path', '_orig_filename', '_tables'):
+            obj.__dict__[attr] = self.__dict__[attr]
+
+        # A new list containing references (rather than a reference to the list)
+        obj.__dict__['_all_nddatas'] = [nd for nd in self._nddata]
+
+        # If we're copying a single-slice AD, then make its NDAstroData a copy
+        # so attributes can be modified without affecting the original
+        if self.is_single:
+            obj.__dict__['_all_nddatas'][self._indices[0]] = copy(self.nddata)
+
+        return obj
+
     def __deepcopy__(self, memo):
         """
         Returns a new instance of this class.

astrodata/utils.py

@@ -266,6 +266,9 @@ def from_shape(value):
     @staticmethod
     def from_string(value):
         """The inverse of __str__, produce a Section object from a string"""
+        # if we were sent None, return None
+        if value is None:
+            return None
         return Section(*[y for x in value.strip("[]").split(",")
                          for start, end in [x.split(":")]
                          for y in (None if start == '' else int(start)-1,

gemini_instruments/f2/adclass.py

@@ -145,11 +145,16 @@ def camera(self, stripID=False, pretty=False):
         # [2022-02-22] The LYOT keyword is now being used to store filters as
         # well as LYOT mask.  Therefore, that keyword cannot be reliably used
         # to return the camera string.
-
-        if self.pixel_scale() > 0.10 and self.pixel_scale() < 0.20:
-            camera = "f/16_G5830"
-        else:
-            camera = self.lyot_stop()
+        camera = self.lyot_stop()
+        if camera is None or not camera.startswith("f/"):
+            focus = self.phu.get("FOCUS", "")
+            if not re.match('^f/\d+$', focus):
+                # Use original pixel scale
+                pixscale = (self.phu["PIXSCALE"] if 'PREPARED' in self.tags
+                            else self.pixel_scale())
+                focus = "f/16" if pixscale > 0.1 else "f/32"
+            # Make up component number (for now) should f/32 happen
+            camera = focus + ("_G5830" if focus == "f/16" else "_G9999")
 
         if camera:
             if stripID or pretty:

gemini_instruments/gemini/adclass.py

@@ -9,6 +9,7 @@
 import math
 import datetime
 import dateutil.parser
+from itertools import chain
 
 import numpy as np
 
@@ -384,17 +385,23 @@ def _dec(self):
 
     def _parse_section(self, keyword, pretty):
         try:
-            value_filter = (str if pretty else Section.from_string)
+            value_filter = lambda x: (x.asIRAFsection() if pretty else x)
             process_fn = lambda x: (None if x is None else value_filter(x))
             # Dummy keyword FULLFRAME returns shape of full data array
             if keyword == 'FULLFRAME':
+                def _encode_shape(shape):
+                    if not shape:
+                        return None
+                    return Section(*chain.from_iterable([(0, npix) for npix in shape[::-1]]))
                 if self.is_single:
-                    sections = '[1:{1},1:{0}]'.format(*self.shape)
+                    sections = _encode_shape(self.shape)
                 else:
-                    sections = ['[1:{1},1:{0}]'.format(*ext.shape)
-                                for ext in self]
+                    sections = [_encode_shape(ext.shape) for ext in self]
             else:
-                sections = self.hdr.get(keyword)
+                if self.is_single:
+                    sections = Section.from_string(self.hdr.get(keyword))
+                else:
+                    sections = [Section.from_string(sec) if sec is not None else None for sec in self.hdr.get(keyword)]
             if self.is_single:
                 return process_fn(sections)
             else:
@@ -2151,9 +2158,11 @@ def empirical_pixel_scale(ext):
 
         pixel_scale_list = []
         for ext in self:
+            if len(ext.shape) < 2:
+                return None
             try:
                 pixel_scale_list.append(3600 * np.sqrt(abs(np.linalg.det(ext.wcs.forward_transform['cd_matrix'].matrix))))
-            except (IndexError, AttributeError):
+            except (IndexError, AttributeError) as iae:
                 scale = empirical_pixel_scale(ext)
                 if scale is not None:
                     pixel_scale_list.append(scale)

gemini_instruments/gnirs/adclass.py

@@ -54,7 +54,7 @@ def _type_spect(self):
             slit = self.phu.get('SLIT', '').lower()
             grat = self.phu.get('GRATING', '')
             prism = self.phu.get('PRISM', '')
-            if slit == 'ifu':
+            if 'ifu' in slit:
                 tags.add('IFU')
             elif ('arcsec' in slit or 'pin' in slit) and 'mm' in grat:
                 if 'MIR' in prism:
@@ -350,6 +350,10 @@ def focal_plane_mask(self, stripID=False, pretty=False):
                 fpm = slit
             elif "XD" in decker:
                 fpm = "{}XD".format(slit)
+            elif "HR-IFU" in slit and "HR-IFU" in decker:
+                fpm = "HR-IFU"
+            elif "LR-IFU" in slit and "LR-IFU" in decker:
+                fpm = "LR-IFU"
             elif "IFU" in slit and "IFU" in decker:
                 fpm = "IFU"
             elif "Acq" in slit and "Acq" in decker:

gemini_instruments/gnirs/tests/test_gnirs.py

@@ -166,5 +166,21 @@ def test_ra_dec_from_text():
     assert ad.target_dec() == pytest.approx(24.345277777777778)
 
 
+# TODO: HR-IFU data doesn't exist yet, add one when available
+EXPECTED_FPMS = [
+    ('N20220918S0040.fits', "LR-IFU"),
+    ("S20061213S0131.fits", "IFU")
+]
+
+
+@pytest.mark.dragons_remote_data
+@pytest.mark.parametrize("filename,expected_fpm", EXPECTED_FPMS)
+def test_ifu_fpm(filename, expected_fpm):
+    path = astrodata.testing.download_from_archive(filename)
+    ad = astrodata.open(path)
+    assert("IFU" in ad.tags)
+    assert(ad.focal_plane_mask(pretty=True) == expected_fpm)
+
+
 if __name__ == "__main__":
     pytest.main()

geminidr/core/parameters_spect.py

@@ -75,23 +75,48 @@ def setDefaults(self):
         del self.grow
 
 
+class createNewApertureConfig(config.Config):
+    aperture = config.Field("Base aperture to offset from", int, None, optional=False)
+    shift = config.Field("Shift (in pixels) to new aperture", float, None, optional=False)
+    aper_upper = config.RangeField("Offset to new upper edge", float, None,
+                                   optional=True, min=0., inclusiveMin=False)
+    aper_lower = config.RangeField("Offset to new lower edge", float, None,
+                                   optional=True, max=0., inclusiveMax=False)
+    suffix = config.Field("Filename suffix", str, "_newApertureCreated", optional=True)
+
+    def validate(self):
+        config.Config.validate(self)
+        if (self.aper_lower and self.aper_upper) or\
+           (not self.aper_lower and not self.aper_upper):
+            pass
+        else:
+            raise ValueError("Both aper_lower and aper_upper must either be "
+                             "specified, or left as None.")
+
+
 class determineDistortionConfig(config.Config):
     suffix = config.Field("Filename suffix", str, "_distortionDetermined", optional=True)
     spatial_order = config.RangeField("Fitting order in spatial direction", int, 3, min=1)
-    spectral_order = config.RangeField("Fitting order in spectral direction", int, 4, min=1)
+    spectral_order = config.RangeField("Fitting order in spectral direction", int, 4, min=0)
     id_only = config.Field("Use only lines identified for wavelength calibration?", bool, False)
     min_snr = config.RangeField("Minimum SNR for peak detection", float, 5., min=3.)
     fwidth = config.RangeField("Feature width in pixels if reidentifying",
-                               float, None, min=2., optional=True)
+                               float, None, min=1., optional=True)
     nsum = config.RangeField("Number of lines to sum", int, 10, min=1)
     step = config.RangeField("Step in rows/columns for tracing", int, 10, min=1)
     max_shift = config.RangeField("Maximum shift per pixel in line position",
                                   float, 0.05, min=0.001, max=0.1)
     max_missed = config.RangeField("Maximum number of steps to miss before a line is lost", int, 5, min=0)
+    min_line_length = config.RangeField("Exclude line traces shorter than this fraction of spatial dimension",
+                                        float, 0., min=0., max=1.)
+    debug_reject_bad = config.Field("Reject lines with suspiciously high SNR (e.g. bad columns)?", bool, True)
     debug = config.Field("Display line traces on image display?", bool, False)
 
+
 class determineSlitEdgesConfig(config.Config):
     suffix = config.Field("Filename suffix", str, "_slitEdgesDetermined", optional=True)
+    spectral_order = config.RangeField("Fitting order in spectral direction",
+                                       int, 3, min=1)
     edges1 = config.ListField("List of left edges of illuminated region(s)",
                               float, default=None, minLength=1,
                               optional=True, single=True)
@@ -101,6 +126,7 @@ class determineSlitEdgesConfig(config.Config):
     debug = config.Field("Plot fits of edges and print extra information?",
                          bool, False)
 
+
 class determineWavelengthSolutionConfig(config.core_1Dfitting_config):
     suffix = config.Field("Filename suffix", str, "_wavelengthSolutionDetermined", optional=True)
     center = config.RangeField("Central row/column to extract", int, None, min=1, optional=True)
@@ -120,8 +146,6 @@ class determineWavelengthSolutionConfig(config.core_1Dfitting_config):
     linelist = config.Field("Filename of arc line list", str, None, optional=True)
     in_vacuo = config.Field("Use vacuum wavelength scale (rather than air)?", bool, False)
     absorption = config.Field("Is feature type absorption?", bool, False)
-    order = config.RangeField("Order of fitting function", int, None, min=0,
-                              optional=True)
     debug_min_lines = config.Field("Minimum number of lines to fit each segment", (str, int), '15,20',
                                    check=list_of_ints_check)
     debug_alternative_centers = config.Field("Try alternative wavelength centers?", bool, False)
@@ -132,6 +156,7 @@ def setDefaults(self):
         del self.grow
         self.niter = 3
 
+
 class distortionCorrectConfig(parameters_generic.calRequirementConfig):
     suffix = config.Field("Filename suffix", str, "_distortionCorrected", optional=True)
     order = config.RangeField("Interpolation order", int, 3, min=0, max=5, inclusiveMax=True)
@@ -217,22 +242,18 @@ class flagCosmicRaysConfig(config.Config):
     spectral_order = config.Field(
         doc="Order for fitting and subtracting object continuum and sky line "
         "models, prior to running the main cosmic ray detection algorithm. "
-        "When None, defaults are used, according to the image size (as in "
-        "the IRAF task gemcrspec). To control which fits are performed, use "
-        "the bkgmodel parameter.",
+        "To control which fits are performed, use the bkgmodel parameter.",
         dtype=int,
         optional=True,
-        default=None,
+        default=9,
     )
     spatial_order = config.Field(
         doc="Order for fitting and subtracting object continuum and sky line "
         "models, prior to running the main cosmic ray detection algorithm. "
-        "When None, defaults are used, according to the image size (as in "
-        "the IRAF task gemcrspec). To control which fits are performed, use "
-        "the bkgmodel parameter.",
+        "To control which fits are performed, use the bkgmodel parameter.",
         dtype=int,
         optional=True,
-        default=None,
+        default=5,
     )
     bkgmodel = config.ChoiceField(
         doc="Set which background model(s) to use, between 'object', "
@@ -420,7 +441,7 @@ class normalizeFlatConfig(config.core_1Dfitting_config):
     center = config.RangeField("Central row/column to extract", int, None, min=1, optional=True)
     nsum = config.RangeField("Number of lines to sum", int, 10, min=1)
     threshold = config.RangeField("Threshold for flagging unilluminated pixels",
-                                  float, 0.1, min=0.01, max=1.0)
+                                  float, 0.01, min=0.01, max=1.0)
     interactive = config.Field("Interactive fitting?", bool, False)
 
     def setDefaults(self):

geminidr/core/primitives_photometry.py

@@ -91,15 +91,26 @@ def addReferenceCatalog(self, adinputs=None, **params):
             try:
                 user_refcat = Table.read(source, format=format)
                 log.stdinfo(f"Successfully read reference catalog {source}")
+            except:
+                log.warning(f"File {source} exists but cannot be read - continuing")
+                return adinputs
+            else:
+                colnames = user_refcat.colnames
                 # Allow raw SExtractor catalog to be used
-                if not {'RAJ2000', 'DEJ2000'}.issubset(user_refcat.colnames):
-                    if {'X_WORLD', 'Y_WORLD'}.issubset(user_refcat.colnames):
+                if not {'RAJ2000', 'DEJ2000'}.issubset(colnames):
+                    if {'X_WORLD', 'Y_WORLD'}.issubset(colnames):
                         log.stdinfo("Renaming X_WORLD, Y_WORLD to RAJ2000, DEJ2000")
                         user_refcat.rename_column('X_WORLD', 'RAJ2000')
                         user_refcat.rename_column('Y_WORLD', 'DEJ2000')
-            except:
-                log.warning(f"File {source} exists but cannot be read - continuing")
-                return adinputs
+                    elif {'RA', 'DEC'}.issubset(colnames):
+                        log.stdinfo("Renaming RA, DEC to RAJ2000, DEJ2000")
+                        user_refcat.rename_column('RA', 'RAJ2000')
+                        user_refcat.rename_column('DEC', 'DEJ2000')
+                if 'Id' in colnames and 'Cat_Id' not in colnames:
+                    log.fullinfo("Renaming Id to Cat_Id")
+                    user_refcat.rename_column('Id', 'Cat_Id')
+                user_refcat.add_column(np.arange(len(user_refcat), dtype=np.int32) + 1,
+                                       name='Id', index=0)
         else:
             user_refcat = None
 

geminidr/core/primitives_register.py

@@ -3,12 +3,9 @@
 #
 #                                                          primitives_register.py
 # ------------------------------------------------------------------------------
-import re
 import math
 import numpy as np
 from astropy.modeling import models
-from astropy import units as u
-from astropy.coordinates import Angle
 
 from copy import deepcopy
 
@@ -247,14 +244,15 @@ def adjustWCSToReference(self, adinputs=None, **params):
             if not use_wcs:
                 ad[0].wcs = deepcopy(adref[0].wcs)
             try:
-                ad[0].wcs.insert_transform(ad[0].wcs.input_frame, transform, after=True)
+                ad[0].wcs.insert_transform(ad[0].wcs.input_frame,
+                                           am.make_serializable(transform), after=True)
             except AttributeError:  # no WCS
-                ad[0].wcs = gWCS([(cf.Frame2D(name="pixels"), transform),
+                ad[0].wcs = gWCS([(cf.Frame2D(name="pixels"), am.make_serializable(transform)),
                                   (cf.Frame2D(name="world"), None)])
 
             # Update X_WORLD and Y_WORLD (RA and DEC) in OBJCAT
             try:
-                x, y = ad[0].OBJCAT['X_IMAGE']-1, ad[0].OBJCAT['Y_IMAGE']-1
+                x, y = ad[0].OBJCAT['X_IMAGE'].data-1, ad[0].OBJCAT['Y_IMAGE'].data-1
             except (AttributeError, KeyError):
                 pass
             else:
@@ -443,9 +441,10 @@ def determineAstrometricSolution(self, adinputs=None, **params):
 
                 if num_matched > 0:
                     # Update OBJCAT (X_WORLD, Y_WORLD)
-                    crpix = (0, 0)  # doesn't matter since we're only doing a shift
+                    crpix = (0, 0)  # doesn't matter since we're only reporting a shift
                     ra0, dec0 = ext.wcs(*crpix)
-                    ext.wcs.insert_transform(ext.wcs.input_frame, transform, after=True)
+                    ext.wcs.insert_transform(ext.wcs.input_frame,
+                                             am.make_serializable(transform), after=True)
                     objcat['X_WORLD'], objcat['Y_WORLD'] = ext.wcs(objcat['X_IMAGE']-1,
                                                                    objcat['Y_IMAGE']-1)