Merge branch 'refs/heads/master' into fluxcal_ee

# Conflicts:
#	python/lvmdrp/main.py

.gitignore

@@ -45,8 +45,11 @@ test_redux-*
 tests/htmlcov
 tests/*coverage*
 *.fits
+*.fits.gz
+*.reg
 *.csv
-
+*.dylib
+*.so
 # Mac OSX
 .DS_Store
 
@@ -60,3 +63,5 @@ de421.bsp
 
 # notebooks
 *checkpoint*
+
+gaia_cache/

README.md

@@ -63,10 +63,21 @@ To install the DRP along with its dependencies, you need to run the following st
     cd lvmdrp
     ```
 
+3. Optional, but recommended: compile c-code for better performance
 
-3. Install the DRP package in the current python environment (see [contributing](#contributing-to-lvm-drp-development) section below for a replacement of this step):
+    ```bash
+    pushd lvmdrp/python/cextern/fast_median/src
+    make
+    popd
+    ```
+    This should have created a fast_median.so (Linux) or fast_median.dylib (MacOS) in the src directory.
+    Set environment variable `LVMDRP_LIB_DIR` if you move the shared library. If you leave it in src/, no need to set this.
+    We will automate this step eventually ;-)
+
+4. Install the DRP package in the current python environment (see [contributing](#contributing-to-lvm-drp-development) section below for a replacement of this step). Make sure you're back in the lvmdrp directory.
 
     ```bash
+    cd lvmdrp
     pip install .
     ```
 

bin/drp

@@ -215,12 +215,11 @@ metacli.add_command(regen)
 cli.add_command(metacli)
 
 
-@cli.command('nightly-cals', short_help='Run the nightly calibration sequence', context_settings=dict(show_default=True))
+@cli.command('calibrations', short_help='Run calibration sequences reduction', context_settings=dict(show_default=True))
 @click.option('-m', '--mjd', type=int, help='an MJD to reduce')
 @click.option('-l', '--mjd-list', type=int, multiple=True, help='a list of specific MJDs to reduce')
 @click.option('-r', '--mjd-range', type=str, help='a range of MJDs to reduce')
 @click.option('-cr', '--reject-cr', is_flag=True, default=False, help='flag to run cosmic rays rejection')
-@click.option('-n', '--use-fiducial-cals', is_flag=True, default=False, help='flag to use fiducial calibrations instead of nightly ones')
 @click.option('-b', '--skip-bias', is_flag=True, default=False, help='skip bias calibrations')
 @click.option('-t', '--skip-trace', is_flag=True, default=False, help='skip fiber traces')
 @click.option('-w', '--skip-wavelength', is_flag=True, default=False, help='skip wavelength calibrations')
@@ -230,62 +229,58 @@ cli.add_command(metacli)
 @click.option('--qrtz-threshold', type=float, default=COUNTS_THRESHOLDS["quartz"], help='count thredhols for tracing with Quartz lamp exposures')
 @click.option('--cent-ncolumns', type=int, default=140, help='number of columns to trace fiber centroids')
 @click.option('--full-ncolumns', type=int, default=40, help='number of columns to full fiber tracing')
+@click.option('-x', '--extract-md', is_flag=True, default=False, help='flag to extract metadata or use cached metadata if exist')
 @click.option('-s', '--skip-done', is_flag=True, default=False, help='flag to skip reduction steps that have already been done')
 @click.option('-c', '--clean-ancillary', is_flag=True, default=False, help='Remove ancillary paths after run')
+@click.option('--fflats-from-to', type=int, default=[None, None], nargs=2, help='Move twilight fiberflats from one MJD to another MJD')
+@click.option('--nightly', is_flag=True, default=False, help='flag to reduce nightly calibration sequence (defaults to long-term)')
 @cloup.constraint(mutually_exclusive, ['mjd', 'mjd_list', 'mjd_range'])
 @cloup.constraint(RequireExactly(1), ['mjd', 'mjd_list', 'mjd_range'])
-def nightly_calibrations(mjd, mjd_list, mjd_range, reject_cr, use_fiducial_cals,
-                         skip_bias, skip_trace, skip_wavelength,
-                         skip_dome, skip_twilight,
-                         ldls_threshold, qrtz_threshold,
-                         cent_ncolumns, full_ncolumns,
-                         skip_done, clean_ancillary):
-    """ Run the nightly calibration sequence """
-
-    # NOTE: skip dome fiberflats while I fix the tracing of the amplitude
-    skip_dome = True
+def calibrations(mjd, mjd_list, mjd_range, reject_cr,
+                 skip_bias, skip_trace, skip_wavelength,
+                 skip_dome, skip_twilight,
+                 ldls_threshold, qrtz_threshold,
+                 cent_ncolumns, full_ncolumns,
+                 extract_md, skip_done, clean_ancillary,
+                 fflats_from_to, nightly):
+    """ Run calibration sequence reduction """
 
     only_cals = {flavor for flavor, skip in zip(["bias", "trace", "wave", "dome", "twilight"],
                                                 [skip_bias, skip_trace, skip_wavelength, skip_dome, skip_twilight]) if not skip}
 
+    # get MJDs for twilight fiberflats copies
+    fflats_from, fflats_to = fflats_from_to
+    if any(fflats_from_to) is None:
+        fflats_from = fflats_to = None
+
     # parse MJDs or
     mjds = parse_mjds(mjd=mjd or mjd_list or mjd_range or None)
     if not isinstance(mjds, (list, tuple)):
         mjds = [mjds]
-    for mjd in mjds:
-        reduce_nightly_sequence(mjd=mjd, reject_cr=reject_cr,
-                                use_fiducial_cals=use_fiducial_cals,
-                                only_cals=only_cals,
-                                counts_thresholds={"ldls": ldls_threshold, "quartz": qrtz_threshold},
-                                cent_guess_ncolumns=cent_ncolumns, trace_full_ncolumns=full_ncolumns,
-                                skip_done=skip_done, keep_ancillary=not clean_ancillary)
-
-# register nightly calibration sequence
-cli.add_command(nightly_calibrations)
 
+    if nightly:
+        for mjd in mjds:
+            reduce_nightly_sequence(mjd=mjd, reject_cr=reject_cr,
+                                    use_longterm_cals=False,
+                                    only_cals=only_cals,
+                                    counts_thresholds={"ldls": ldls_threshold, "quartz": qrtz_threshold},
+                                    cent_guess_ncolumns=cent_ncolumns, trace_full_ncolumns=full_ncolumns,
+                                    extract_metadata=extract_md, skip_done=skip_done,
+                                    keep_ancillary=not clean_ancillary,
+                                    fflats_from=fflats_from if fflats_to == mjd else None)
+    else:
+        for mjd in mjds:
+            reduce_longterm_sequence(mjd=mjd, reject_cr=reject_cr,
+                                    use_longterm_cals=True,
+                                    only_cals=only_cals,
+                                    counts_thresholds={"ldls": ldls_threshold, "quartz": qrtz_threshold},
+                                    cent_guess_ncolumns=cent_ncolumns, trace_full_ncolumns=full_ncolumns,
+                                    extract_metadata=extract_md,
+                                    skip_done=skip_done, keep_ancillary=not clean_ancillary,
+                                    fflats_from=fflats_from if fflats_to == mjd else None)
 
-@cli.command('long-term-cals', short_help='Run the long-term calibration sequence')
-@click.option('-m', '--mjd', type=int, required=True, help='the MJD to reduce')
-@click.option('-cr', '--reject-cr', is_flag=True, default=False, help='flag to run cosmic rays rejection')
-@click.option('-b', '--skip-bias', is_flag=True, default=False, help='skip bias calibrations')
-@click.option('-t', '--skip-trace', is_flag=True, default=False, help='skip fiber traces')
-@click.option('-w', '--skip-wavelength', is_flag=True, default=False, help='skip wavelength calibrations')
-@click.option('-f', '--skip-fiberflat', is_flag=True, default=False, help='skip fiberflats')
-@click.option('-s', '--skip-done', is_flag=True, default=False, help='flag to skip reduction steps that have already been done')
-@click.option('-c', '--clean-ancillary', is_flag=True, default=False, help='flag to keep ancillary files')
-def long_term_calibrations(mjd, reject_cr, skip_bias, skip_trace, skip_wavelength,
-                           skip_fiberflat, skip_done, clean_ancillary):
-    """ Run the long-term calibration sequence """
-
-    only_cals = {flavor for flavor, skip in zip(["bias", "trace", "wave", "twilight"],
-                                                [skip_bias, skip_trace, skip_wavelength, skip_fiberflat]) if not skip}
-    reduce_longterm_sequence(mjd=mjd, reject_cr=reject_cr,
-                             use_fiducial_cals=True,
-                             only_cals=only_cals,
-                             skip_done=skip_done, keep_ancillary=not clean_ancillary)
-
-# register long-term calibration sequence
-cli.add_command(long_term_calibrations)
+# register calibration sequence
+cli.add_command(calibrations)
 
 
 @cli.command('fix-pixel-shifts', short_help='Fix the electronic pixel shifts in the raw data')
@@ -364,8 +359,10 @@ def get_calibs(kind, mjd, camera):
 @click.option('-w', '--walltime', type=str, default="24:00:00", help='the time for which the job is allowed to run')
 @click.option('-a', '--alloc', type=click.Choice(['sdss-np', 'sdss-kp']), default='sdss-np', help='which partition to use')
 @click.option('-s', '--submit', is_flag=True, type=bool, default=True, help='flag to submit the job or not')
-@click.option('--run-calibs', is_flag=True, type=bool, default=False, help="run long-term calibrations only")
-def cluster(mjd_list, mjd_range, expnum, exp_list, exp_range, exp_file, nodes, ppn, walltime, alloc, submit, run_calibs):
+@click.option('--run-calibs', is_flag=True, type=bool, default=False, help="run calibration sequences only")
+@click.option('--drp-options', type=str, default="", help="pass options to drp run command")
+def cluster(mjd_list, mjd_range, expnum, exp_list, exp_range, exp_file, nodes, ppn, walltime,
+            alloc, submit, run_calibs, drp_options):
     """ Submit a Utah cluster job to batch run the DRP by MJD """
 
     # filter the mjds
@@ -375,7 +372,8 @@ def cluster(mjd_list, mjd_range, expnum, exp_list, exp_range, exp_file, nodes, p
         mjds = parse_mjds(mjds)
 
     # submit the cluster job
-    run_cluster(mjds=mjds, expnums=expnums, nodes=nodes, ppn=ppn, walltime=walltime, alloc=alloc, submit=submit, run_calibs=run_calibs)
+    run_cluster(mjds=mjds, expnums=expnums, nodes=nodes, ppn=ppn, walltime=walltime, alloc=alloc, submit=submit,
+                run_calibs=run_calibs, drp_options=drp_options)
 
 
 if __name__ == "__main__":

bin/pix2wave

@@ -151,6 +151,11 @@ class SolveWavelength:
         # read given image
         log.info(f"loading extracted RSS '{self.extracted_fits}'")
         spectra = rss.loadRSS(self.extracted_fits)
+        select = spectra._mask|np.isnan(spectra._data)
+        spectra._data[select] = 0.0
+        spectra._error[select] = np.inf
+        log.info("subtracting continuum")
+        spectra, _, _ = spectra.subtract_continuum(niter=3)
         # extract relevant metadata
         self.lamps = [lamp.lower() for lamp in ARC_LAMPS if spectra._header.get(lamp, "OFF") == "ON"]
         self.camera = spectra._header["CCD"]
@@ -225,11 +230,13 @@ class SolveWavelength:
         (self.line,) = self.ax_pix.step(
             self.xdata,
             self.ydata,
+            where="mid",
             lw=1,
             picker=True,
             pickradius=self.aperture,
             color="k",
         )
+        self.ax_pix.set_yscale("log")
 
         # fix table with cross-matching
         if self.table_exist and self.correct_lines:
@@ -246,7 +253,7 @@ class SolveWavelength:
         # create guess wavelength function
         poly = self._fit_wavelength()
         waves_guess = poly(self.xdata)
-        self.ax_wav.step(waves_guess, self.ydata, lw=1, color="k")
+        self.ax_wav.step(waves_guess, self.ydata, where="mid", lw=1, color="k")
         for _, ref_wave, name in self.reflines:
             self.ax_wav.axvline(ref_wave, lw=1, color="tab:red")
             self.ax_wav.text(
@@ -265,6 +272,7 @@ class SolveWavelength:
         self.ax_wav.set_xlim(*self.wave_range)
         self.ax_wav.set_xlabel("wavelength (A)")
         self.ax_wav.set_ylabel("counts (e-/pix)")
+        self.ax_wav.set_yscale("log")
 
         # adding known peaks
         log.info(f"displaying {self.pixmap.shape[0]} known lines (including masked)")