Merge pull request #1710 from pypeit/mdm_r4k

Adds R4K support for MDM OSMOS

doc/releases/1.14.1dev.rst

@@ -24,6 +24,7 @@ Functionality/Performance Improvements and Additions
 - Add a sensible error message to the pypeit Spectrum1D loaders in the event a
   user inadvertently tries to use Spectrum1D instead of SpectrumList for a
   ``spec1d`` file.
+- Add support for the R4K detector for MDM OSMOS
 
 Instrument-specific Updates
 ---------------------------

doc/spectrographs/mdm_osmos.rst

@@ -0,0 +1,20 @@
+*********
+MDM OSMOS
+*********
+
+
+Overview
+========
+
+This file summarizes items related to MDM OSMOS
+
+
+R4K
+===
+
+It is common for only a portion of 
+the data one receives to have been
+windowed.  You may therefore need to 
+window the rest.  There is a Notebook
+in the DevSuite that shows an example of
+how to do this.

pypeit/core/procimg.py

@@ -16,7 +16,6 @@
 
 from pypeit import msgs
 from pypeit import utils
-from pypeit.core import parse
 
 
 # NOTE: This is slower than utils.rebin_evlist by a factor of ~2, but avoids an
@@ -589,7 +588,14 @@ def rect_slice_with_mask(image, mask, mask_val=1):
 def subtract_overscan(rawframe, datasec_img, oscansec_img, method='savgol', params=[5,65],
                       var=None):
     """
-    Subtract overscan.
+    Subtract the overscan
+
+    Possible values of ``method``:
+        - polynomial: Fit a polynomial to the overscan region and subtract it.
+        - savgol: Use a Savitzky-Golay filter to fit the overscan region and
+            subtract it.
+        - median: Use the median of the overscan region to subtract it.
+        - odd_even: Use the median of the odd and even rows/columns to subtract (MDM/OSMOS)
 
     Args:
         rawframe (`numpy.ndarray`_):
@@ -622,11 +628,11 @@ def subtract_overscan(rawframe, datasec_img, oscansec_img, method='savgol', para
     Returns:
         :obj:`tuple`: The input frame with the overscan region subtracted and an
         estimate of the variance in the overscan subtraction; both have the same
-        shape as the input ``rawframe``.  If ``var`` is no provided, the 2nd
+        shape as the input ``rawframe``.  If ``var`` is not provided, the 2nd
         returned object is None.
     """
     # Check input
-    if method.lower() not in ['polynomial', 'savgol', 'median']:
+    if method.lower() not in ['polynomial', 'savgol', 'median', 'odd_even']:
         msgs.error(f'Unrecognized overscan subtraction method: {method}')
     if rawframe.ndim != 2:
         msgs.error('Input raw frame must be 2D.')
@@ -671,6 +677,7 @@ def subtract_overscan(rawframe, datasec_img, oscansec_img, method='savgol', para
             # to the error in the mean
             osvar = np.pi/2*(np.sum(osvar)/osvar.size**2 if method.lower() == 'median' 
                              else np.sum(osvar, axis=compress_axis)/osvar.shape[compress_axis]**2)
+        # Method time
         if method.lower() == 'polynomial':
             # TODO: Use np.polynomial.polynomial.polyfit instead?
             c = np.polyfit(np.arange(osfit.size), osfit, params[0])
@@ -683,12 +690,25 @@ def subtract_overscan(rawframe, datasec_img, oscansec_img, method='savgol', para
             if var is not None:
                 _var[data_slice] = osvar
             continue
+        elif method.lower() == 'odd_even':
+            ossub = np.zeros_like(osfit)
+            # Odd/even
+            if compress_axis == 1:
+                odd = np.median(overscan[:,1::2], axis=compress_axis)
+                even = np.median(overscan[:,0::2], axis=compress_axis)
+                # Do it
+                no_overscan[data_slice][:,1::2] -= odd[:,None]
+                no_overscan[data_slice][:,0::2] -= even[:,None]
+            else:
+                msgs.error('Not ready for this approach, please contact the Developers')
+
 
         # Subtract along the appropriate axis
         no_overscan[data_slice] -= (ossub[:, None] if compress_axis == 1 else ossub[None, :])
         if var is not None:
             _var[data_slice] = (osvar[:,None] if compress_axis == 1 else osvar[None,:])
 
+    # Return
     return no_overscan, _var
 
 

pypeit/core/wavecal/autoid.py

@@ -996,7 +996,6 @@ def set_fwhm(par, measured_fwhm=None, verbose=False):
     return fwhm
 
 
-# TODO: Docstring is missing many arguments
 def full_template(spec, lamps, par, ok_mask, det, binspectral, nsnippet=2, slit_ids=None,
                   measured_fwhms=None, debug_xcorr=False, debug_reid=False,
                   x_percentile=50., template_dict=None, debug=False, 
@@ -1038,6 +1037,16 @@ def full_template(spec, lamps, par, ok_mask, det, binspectral, nsnippet=2, slit_
         Passed to reidentify to reduce the dynamic range of arc line amplitudes
     template_dict : dict, optional
         Dict containing tempmlate items, largely for development
+    nonlinear_counts : float, optional
+        For arc line detection: Arc lines above this saturation threshold
+        are not used in wavelength solution fits because they cannot be
+        accurately centroided. Defaults to 1e10.
+    debug : bool, optional
+        Show plots useful for debugging
+    debug_xcorr : bool, optional
+        Show plots useful for debugging the cross-correlation
+    debug_reid : bool, optional
+        Show plots useful for debugging the reidentification
 
     Returns
     -------

pypeit/images/rawimage.py

@@ -1017,6 +1017,8 @@ def subtract_overscan(self, force=False):
         """
         Analyze and subtract the overscan from the image
 
+        If this is a mosaic, loop over the individual detectors
+
         Args:
             force (:obj:`bool`, optional):
                 Force the image to be overscan subtracted, even if the step log
@@ -1036,12 +1038,11 @@ def subtract_overscan(self, force=False):
         for i in range(self.nimg):
             # Subtract the overscan.  var is the variance in the overscan
             # subtraction.
-            _os_img[i], var[i] = procimg.subtract_overscan(self.image[i], self.datasec_img[i],
-                                                           self.oscansec_img[i],
-                                                           method=self.par['overscan_method'],
-                                                           params=self.par['overscan_par'],
-                                                           var=None if self.rn2img is None 
-                                                                else self.rn2img[i])
+            _os_img[i], var[i] = procimg.subtract_overscan(
+                self.image[i], self.datasec_img[i], self.oscansec_img[i],
+                method=self.par['overscan_method'],
+                params=self.par['overscan_par'],
+                var=None if self.rn2img is None else self.rn2img[i])
         self.image = np.array(_os_img)
         # Parse the returned value
         if self.rn2img is not None:

pypeit/par/pypeitpar.py

@@ -459,7 +459,7 @@ def valid_overscan_methods():
         """
         Return the valid overscan methods.
         """
-        return ['polynomial', 'savgol', 'median']
+        return ['polynomial', 'savgol', 'median', 'odd_even']
 
     @staticmethod
     def valid_combine_methods():

pypeit/spectrographs/mdm_osmos.py

@@ -80,6 +80,7 @@ def default_pypeit_par(cls):
         """
         par = super().default_pypeit_par()
 
+
         # Ignore PCA
         par['calibrations']['slitedges']['sync_predict'] = 'nearest'
 
@@ -233,3 +234,135 @@ def check_frame_type(self, ftype, fitstbl, exprng=None):
         msgs.warn('Cannot determine if frames are of type {0}.'.format(ftype))
         return np.zeros(len(fitstbl), dtype=bool)
 
+
+class MDMOSMOSR4KSpectrograph(MDMOSMOSMDM4KSpectrograph):
+    """
+    Child to handle MDM OSMOS R4K instrument+detector
+    """
+    ndet = 1
+    name = 'mdm_osmos_r4k'
+    camera = 'R4K'
+    url = 'https://www.astronomy.ohio-state.edu/martini.10/osmos/'
+    header_name = 'OSMOS'
+    supported = True
+    comment = 'MDM OSMOS spectrometer for the red. Requires calibrations windowed down to the science frame.'
+
+    def get_detector_par(self, det, hdu=None):
+        """
+        Return metadata for the selected detector.
+
+        THIS IS FOR WINDOWED SCIENCE FRAMES
+        AND WE ARE HACKING THE CALIBS
+
+        Args:
+            det (:obj:`int`):
+                1-indexed detector number.
+            hdu (`astropy.io.fits.HDUList`_, optional):
+                The open fits file with the raw image of interest.  If not
+                provided, frame-dependent parameters are set to a default.
+
+        Returns:
+            :class:`~pypeit.images.detector_container.DetectorContainer`:
+            Object with the detector metadata.
+        """
+        # Detector 1
+        detector_dict = dict(
+            binning         = '1,1' if hdu is None 
+                                    else self.get_meta_value(self.get_headarr(hdu), 'binning'),
+            det=1,
+            dataext         = 0,
+            specaxis        = 1,
+            specflip        = True,
+            spatflip        = False,
+            xgap            = 0.,
+            ygap            = 0.,
+            ysize           = 1.,
+            platescale      = 0.273,
+            mincounts       = -1e10,
+            darkcurr        = 0.0,
+            saturation      = 65535.,
+            nonlinear       = 0.86,
+            numamplifiers   = 4,
+            gain            = np.atleast_1d([2.2, 2.2, 2.2, 2.2]),
+            ronoise         = np.atleast_1d([3.0, 3.0, 3.0, 3.0]),
+            datasec         = np.atleast_1d(['[:524,33:2064]', '[524:,33:2064]',
+                '[:524, 2065:4092', '[524:, 2065:4092']),
+            oscansec        = np.atleast_1d(['[:524, 1:32]', '[524:, 1:32]',
+                '[:524, 4129:]', '[524:, 4129:]']),
+        )
+        # Return
+        return detector_container.DetectorContainer(**detector_dict)
+
+    def check_frame_type(self, ftype, fitstbl, exprng=None):
+        """
+        Check for frames of the provided type.
+
+        Args:
+            ftype (:obj:`str`):
+                Type of frame to check. Must be a valid frame type; see
+                frame-type :ref:`frame_type_defs`.
+            fitstbl (`astropy.table.Table`_):
+                The table with the metadata for one or more frames to check.
+            exprng (:obj:`list`, optional):
+                Range in the allowed exposure time for a frame of type
+                ``ftype``. See
+                :func:`pypeit.core.framematch.check_frame_exptime`.
+
+        Returns:
+            `numpy.ndarray`_: Boolean array with the flags selecting the
+            exposures in ``fitstbl`` that are ``ftype`` type frames.
+        """
+        good_exp = framematch.check_frame_exptime(fitstbl['exptime'], exprng)
+        if ftype in ['science', 'standard']:
+            return good_exp & (fitstbl['idname'] == 'OBJECT')
+        if ftype == 'bias':
+            return good_exp & (fitstbl['idname'] == 'zero')
+        if ftype in ['pixelflat', 'trace']:
+            return good_exp & (fitstbl['lampstat01'] == 'Flat') & (fitstbl['idname'] == 'FLAT')
+        if ftype in ['pinhole', 'dark']:
+            # Don't type pinhole or dark frames
+            return np.zeros(len(fitstbl), dtype=bool)
+        if ftype in ['arc','tilt']:
+            return good_exp & (fitstbl['idname'] == 'COMP') 
+        msgs.warn('Cannot determine if frames are of type {0}.'.format(ftype))
+        return np.zeros(len(fitstbl), dtype=bool)
+
+    @classmethod
+    def default_pypeit_par(cls):
+        """
+        Return the default parameters to use for this instrument.
+        
+        Returns:
+            :class:`~pypeit.par.pypeitpar.PypeItPar`: Parameters required by
+            all of PypeIt methods.
+        """
+        par = super().default_pypeit_par()
+
+        # Do not require bias frames
+        turn_off = dict(use_biasimage=False, overscan_method='odd_even')
+        par.reset_all_processimages_par(**turn_off)
+
+        # Ignore PCA
+        par['calibrations']['slitedges']['sync_predict'] = 'nearest'
+        # Bound the detector with slit edges if no edges are found
+        par['calibrations']['slitedges']['bound_detector'] = True
+
+        # Set pixel flat combination method
+        par['calibrations']['pixelflatframe']['process']['combine'] = 'median'
+        # Wavelength calibration methods
+        par['calibrations']['wavelengths']['method'] = 'full_template'
+        #par['calibrations']['wavelengths']['method'] = 'reidentify'
+        par['calibrations']['wavelengths']['lamps'] = ['HgI', 'NeI']
+        par['calibrations']['wavelengths']['reid_arxiv'] = 'mdm_osmos_r4k.fits'
+        par['calibrations']['wavelengths']['sigdetect'] = 5.0
+        par['calibrations']['wavelengths']['nsnippet'] = 1
+        par['calibrations']['wavelengths']['fwhm_fromlines'] = True
+        # Set the default exposure time ranges for the frame typing
+        par['calibrations']['biasframe']['exprng'] = [None, 1]
+        par['calibrations']['darkframe']['exprng'] = [999999, None]     # No dark frames
+        par['calibrations']['pinholeframe']['exprng'] = [999999, None]  # No pinhole frames
+        par['calibrations']['arcframe']['exprng'] = [None, None]  # Long arc exposures on this telescope
+        par['calibrations']['standardframe']['exprng'] = [None, 120]
+        par['scienceframe']['exprng'] = [90, None]
+
+        return par

pypeit/wavecalib.py

@@ -658,7 +658,7 @@ def build_wv_calib(self, arccen, method, skip_QA=False,
                                              self.binspectral, slit_ids=self.slits.slitord_id,
                                              measured_fwhms=self.measured_fwhms,
                                              nonlinear_counts=self.nonlinear_counts,
-                                             nsnippet=self.par['nsnippet'])
+                                             nsnippet=self.par['nsnippet'])#,
                                              #debug=True, debug_reid=True, debug_xcorr=True)
         elif self.par['method'] == 'echelle':
             # Echelle calibration files