Merge pull request #82 from gbrammer/extend-wavelength-range

Extend NIRSpec wavelength ranges

.flake8

@@ -0,0 +1,7 @@
+[flake8]
+count = true
+select = E9,F63,F7,F82
+show-source = true
+statistics = true
+exclude = .git,build,docs,scripts,binder
+

.github/workflows/python-package-with-numba.yml

@@ -48,7 +48,7 @@ jobs:
         pip install flake8
     - name: Lint with flake8
       run: |
-        flake8 . --count --select=E9,F63,F7,F82 --show-source --statistics
+        flake8 .
   test:
     runs-on: ${{ matrix.os }}
     strategy:

.github/workflows/python-package.yml

@@ -48,7 +48,7 @@ jobs:
         pip install flake8
     - name: Lint with flake8
       run: |
-        flake8 . --count --select=E9,F63,F7,F82 --show-source --statistics
+        flake8 .
   test:
     runs-on: ${{ matrix.os }}
     strategy:

docs/msaexp/api.rst

@@ -3,6 +3,9 @@ Additional API
 .. automodapi:: msaexp.pipeline
     :no-inheritance-diagram:
 
+.. automodapi:: msaexp.pipeline_extended
+    :no-inheritance-diagram:
+
 .. automodapi:: msaexp.utils
     :no-inheritance-diagram:
 

msaexp/data/fwhm_data.py

@@ -112,5 +112,36 @@
 5.38 1.859
 """
 
+# Measured again from 1493 with multiple Gaussians
+fwhm_data_2024b = """# wave fwhm_pix
+0.50   1.24
+0.58   1.24
+0.66   1.25
+0.79   1.26
+0.99   1.27
+1.29   1.29
+1.65   1.29
+2.00   1.30
+2.32   1.31
+2.61   1.31
+2.87   1.33
+3.11   1.33
+3.34   1.35
+3.55   1.37
+3.75   1.39
+3.94   1.43
+4.12   1.46
+4.29   1.51
+4.46   1.55
+4.62   1.60
+4.77   1.64
+4.92   1.67
+5.06   1.70
+5.20   1.76
+5.34   1.82
+5.60   1.90
+"""
+
 # Set default
-fwhm_data = fwhm_data_update
+# fwhm_data = fwhm_data_update
+fwhm_data = fwhm_data_2024b

msaexp/drizzle.py

@@ -908,8 +908,8 @@ def get_xlimits_from_lines(
     """
     Find emission lines in 2D spectrum
 
-    Parameters:
-    -----------
+    Parameters
+    ----------
     hdul : HDUList
         The HDUList object containing the 2D spectrum.
 
@@ -933,8 +933,8 @@ def get_xlimits_from_lines(
     verbose : bool, optional
         Whether to print verbose output. Default is True.
 
-    Returns:
-    --------
+    Returns
+    -------
     xlim : tuple
         A tuple containing the x-limits of the emission lines.
 
@@ -1467,7 +1467,12 @@ def extract_from_hdul(
         )
 
     if prf_center is None:
-        prf_center = sci.header["SRCYPIX"] - (sci.data.shape[0] - 1) / 2.0
+        if "SRCYPIX" in sci.header:
+            y0 = sci.header["SRCYPIX"]
+        else:
+            y0 = (sci.data.shape[0] - 1) / 2
+
+        prf_center = y0 - (sci.data.shape[0] - 1) / 2.0
 
     _data = make_optimal_extraction(
         waves,

msaexp/pipeline.py

@@ -410,6 +410,137 @@ def read_data(self, verbose=True):
             )
 
 
+def exposure_oneoverf(file, fix_rows=False, skip_completed=True, **kwargs):
+    """
+    Remove column-averaged 1/f striping
+
+    Parameters
+    ----------
+    file : str
+        Exposure (rate.fits) filename
+
+    fix_rows : bool
+        Apply 1/f correction to detector rows, as well as columns
+
+    skip_completed : bool
+        Skip steps that have already been completed
+
+    Returns
+    -------
+    status : bool
+        False if ``skip_completed`` and ``ONEFEXP`` keyword found, True if
+        executed without exception.
+
+    """
+    with pyfits.open(file) as im:
+        if "ONEFEXP" in im[0].header:
+            if im[0].header["ONEFEXP"] & skip_completed:
+                return False
+
+    jwst_utils.exposure_oneoverf_correction(
+        file, erode_mask=False, in_place=True, axis=0, deg_pix=256
+    )
+
+    if fix_rows:
+        jwst_utils.exposure_oneoverf_correction(
+            file, erode_mask=False, in_place=True, axis=1, deg_pix=2048
+        )
+
+    return True
+
+
+def exposure_detector_effects(
+    file, fix_rows=False, scale_rnoise=True, skip_completed=True, **kwargs
+):
+    """
+    Remove 1/f striping, bias pedestal offset and rescale RNOISE extension
+
+    Parameters
+    ----------
+    file : str
+        Exposure (rate.fits) filename
+
+    scale_rnoise : bool
+        Calculate the RNOISE scaling
+
+    skip_completed : bool
+        Skip steps that have already been completed
+
+    Returns
+    -------
+    status : bool
+        False if ``skip_completed`` and ``ONEFEXP`` keyword found, True if
+        executed without exception.
+
+    """
+
+    status = exposure_oneoverf(
+        file, fix_rows=fix_rows, skip_completed=skip_completed, **kwargs
+    )
+
+    with pyfits.open(file, mode="update") as im:
+        # bias
+        dq = (im["DQ"].data & 1025) == 0
+
+        if im[0].header["DETECTOR"] == "NRS2":
+            dq[:, :1400] = False
+        else:
+            dq[:, 1400:] = False
+
+        if ("MASKBIAS" in im[0].header) & skip_completed:
+            bias_level = im[0].header["MASKBIAS"]
+            msg = f"msaexp.preprocess : {file}  bias offset ="
+            msg += f" {bias_level:7.3f} (from MASKBIAS)"
+            utils.log_comment(
+                utils.LOGFILE, msg, verbose=True, show_date=False
+            )
+        else:
+            bias_level = np.nanmedian(im["SCI"].data[dq])
+            msg = f"msaexp.preprocess : {file}  bias offset ="
+            msg += f" {bias_level:7.3f}"
+            utils.log_comment(
+                utils.LOGFILE, msg, verbose=True, show_date=False
+            )
+
+            im["SCI"].data -= bias_level
+            im[0].header["MASKBIAS"] = bias_level, "Bias level"
+            im[0].header["MASKNPIX"] = (
+                dq.sum(),
+                "Number of pixels used for bias level",
+            )
+
+        if scale_rnoise:
+
+            if ("SCLREADN" in im[0].header) & skip_completed:
+                rms = im[0].header["SCLREADN"]
+                msg = f"msaexp.preprocess : {file}    rms scale ="
+                msg += f" {rms:>7.2f} (from SCLREADN)"
+                utils.log_comment(
+                    utils.LOGFILE, msg, verbose=True, show_date=False
+                )
+            else:
+                resid = im["SCI"].data / np.sqrt(im["VAR_RNOISE"].data)
+                rms = utils.nmad(resid[dq])
+                msg = f"msaexp.preprocess : {file}    rms scale ="
+                msg += f"{rms:>7.2f}"
+                utils.log_comment(
+                    utils.LOGFILE, msg, verbose=True, show_date=False
+                )
+
+                im[0].header["SCLREADN"] = rms, "RNOISE Scale factor"
+
+                im["VAR_RNOISE"].data *= rms**2
+
+                im[0].header["SCLRNPIX"] = (
+                    dq.sum(),
+                    "Number of pixels used for rnoise scale",
+                )
+
+        im.flush()
+
+    return True
+
+
 class NirspecPipeline:
     def __init__(
         self,
@@ -735,83 +866,15 @@ def preprocess(
             snowball_dilate=24,
         )
 
-        # 1/f correction
+        # 1/f, bias & rnoise
         for file in self.files:
-            with pyfits.open(file) as im:
-                if "ONEFEXP" in im[0].header:
-                    if im[0].header["ONEFEXP"] & skip_completed:
-                        continue
-
-            jwst_utils.exposure_oneoverf_correction(
-                file, erode_mask=False, in_place=True, axis=0, deg_pix=256
+            exposure_detector_effects(
+                file,
+                fix_rows=fix_rows,
+                scale_rnoise=scale_rnoise,
+                skip_completed=skip_completed,
             )
 
-            if fix_rows:
-                jwst_utils.exposure_oneoverf_correction(
-                    file, erode_mask=False, in_place=True, axis=1, deg_pix=2048
-                )
-
-        # bias
-        for file in self.files:
-            with pyfits.open(file, mode="update") as im:
-                dq = (im["DQ"].data & 1025) == 0
-
-                if im[0].header["DETECTOR"] == "NRS2":
-                    dq[:, :1400] = False
-                else:
-                    dq[:, 1400:] = False
-
-                if ("MASKBIAS" in im[0].header) & skip_completed:
-                    bias_level = im[0].header["MASKBIAS"]
-                    msg = f"msaexp.preprocess : {file}  bias offset ="
-                    msg += f" {bias_level:7.3f} (from MASKBIAS)"
-                    utils.log_comment(
-                        utils.LOGFILE, msg, verbose=True, show_date=False
-                    )
-                else:
-                    bias_level = np.nanmedian(im["SCI"].data[dq])
-                    msg = f"msaexp.preprocess : {file}  bias offset ="
-                    msg += f" {bias_level:7.3f}"
-                    utils.log_comment(
-                        utils.LOGFILE, msg, verbose=True, show_date=False
-                    )
-
-                    im["SCI"].data -= bias_level
-                    im[0].header["MASKBIAS"] = bias_level, "Bias level"
-                    im[0].header["MASKNPIX"] = (
-                        dq.sum(),
-                        "Number of pixels used for bias level",
-                    )
-
-                if scale_rnoise:
-
-                    if ("SCLREADN" in im[0].header) & skip_completed:
-                        rms = im[0].header["SCLREADN"]
-                        msg = f"msaexp.preprocess : {file}    rms scale ="
-                        msg += f" {rms:>7.2f} (from SCLREADN)"
-                        utils.log_comment(
-                            utils.LOGFILE, msg, verbose=True, show_date=False
-                        )
-                    else:
-                        resid = im["SCI"].data / np.sqrt(im["VAR_RNOISE"].data)
-                        rms = utils.nmad(resid[dq])
-                        msg = f"msaexp.preprocess : {file}    rms scale ="
-                        msg += f"{rms:>7.2f}"
-                        utils.log_comment(
-                            utils.LOGFILE, msg, verbose=True, show_date=False
-                        )
-
-                        im[0].header["SCLREADN"] = rms, "RNOISE Scale factor"
-
-                        im["VAR_RNOISE"].data *= rms**2
-
-                        im[0].header["SCLRNPIX"] = (
-                            dq.sum(),
-                            "Number of pixels used for rnoise scale",
-                        )
-
-                im.flush()
-
         return True
 
     def run_jwst_pipeline(
@@ -1571,8 +1634,7 @@ def extract_spectrum(
         center2d=False,
         trace_sign=1,
         min_dyoffset=0.2,
-        **kwargs,
-    ):
+        **kwargs):
         """
         Main function for extracting 2D/1D spectra from individual slitlets
 
@@ -1637,8 +1699,8 @@ def extract_spectrum(
         min_dyoffset: float
             The minimum value for the y offset. Default is 0.2.
 
-        Returns:
-        ----------
+        Returns
+        -------
         If 'key' is not found in 'slitlets', returns None.
 
         Else:
@@ -2487,17 +2549,16 @@ def make_summary_tables(root="msaexp", zout=None):
     """
     Make a summary table of all extracted sources
 
-    Parameters:
-    -----------
-
+    Parameters
+    ----------
     root : str
         The root directory where the data is stored. Default is "msaexp".
 
     zout : astropy.table.Table
         Optional table containing photometric redshift information. Default is None.
 
-    Returns:
-    --------
+    Returns
+    -------
     tabs : list
         List of astropy tables containing the extracted source information.