DM-31489: Update StrayLightData to use FitsGenericFormatter with a deferred data set

python/lsst/obs/subaru/strayLight/formatter.py

@@ -19,24 +19,23 @@
 # You should have received a copy of the GNU General Public License
 # along with this program.  If not, see <http://www.gnu.org/licenses/>.
 
-from lsst.daf.butler.formatters.file import FileFormatter
+from lsst.obs.base.formatters.fitsGeneric import FitsGenericFormatter
 
 from .yStrayLight import SubaruStrayLightData
 
 __all__ = ("SubaruStrayLightDataFormatter",)
 
 
-class SubaruStrayLightDataFormatter(FileFormatter):
+class SubaruStrayLightDataFormatter(FitsGenericFormatter):
     """Gen3 Butler Formatters for HSC y-band stray light correction data.
     """
-    extension = ".fits"
 
     extension = ".fits"
 
     def _readFile(self, path, pytype=None):
-        # Docstring inherited from FileFormatter._readFile.
-        return SubaruStrayLightData(path)
+        # Docstring inherited from FitsGenericFormatter._readFile.
+        return SubaruStrayLightData.readFits(path)
 
     def _writeFile(self, inMemoryDataset):
-        # Docstring inherited from FileFormatter._writeFile.
+        # Docstring inherited from FitsGenericFormatter._writeFile.
         raise NotImplementedError("SubaruStrayLightData cannot be written directly.")

python/lsst/obs/subaru/strayLight/yStrayLight.py

@@ -24,6 +24,7 @@
 import scipy.interpolate
 
 from lsst.geom import Angle, degrees
+from lsst.daf.butler import DeferredDatasetHandle
 from lsst.ip.isr.straylight import StrayLightConfig, StrayLightTask, StrayLightData
 
 from . import waveletCompression
@@ -33,7 +34,6 @@
 BAD_THRESHOLD = 500  # Threshold for identifying bad pixels in the reconstructed dark image
 
 
-# TODO DM-16805: This doesn't match the rest of the obs_subaru/ISR code.
 class SubaruStrayLightTask(StrayLightTask):
     """Remove stray light in the y-band
 
@@ -49,6 +49,7 @@ class SubaruStrayLightTask(StrayLightTask):
     This Task retrieves the appropriate dark, uncompresses it and
     uses it to remove the stray light from an exposure.
     """
+
     ConfigClass = StrayLightConfig
 
     def readIsrData(self, dataRef, rawExposure):
@@ -59,7 +60,7 @@ def readIsrData(self, dataRef, rawExposure):
         if not self.check(rawExposure):
             return None
 
-        return SubaruStrayLightData(dataRef.get("yBackground_filename")[0])
+        return SubaruStrayLightData.readFits(dataRef.get("yBackground_filename")[0])
 
     def check(self, exposure):
         # Docstring inherited from StrayLightTask.check.
@@ -101,6 +102,10 @@ def run(self, exposure, strayLightData):
         if strayLightData is None:
             raise RuntimeError("No strayLightData supplied for correction.")
 
+        if isinstance(strayLightData, DeferredDatasetHandle):
+            # Get the deferred object.
+            strayLightData = strayLightData.get()
+
         exposureMetadata = exposure.getMetadata()
         detId = exposure.getDetector().getId()
         if self.config.doRotatorAngleCorrection:
@@ -141,7 +146,7 @@ def run(self, exposure, strayLightData):
 
 
 class SubaruStrayLightData(StrayLightData):
-    """Lazy-load object that reads and integrates the wavelet-compressed
+    """Object that reads and integrates the wavelet-compressed
     HSC y-band stray-light model.
 
     Parameters
@@ -150,8 +155,20 @@ class SubaruStrayLightData(StrayLightData):
         Full path to a FITS files containing the stray-light model.
     """
 
-    def __init__(self, filename):
-        self._filename = filename
+    @classmethod
+    def readFits(cls, filename, **kwargs):
+        calib = cls()
+        hdulist = fits.open(filename)
+
+        calib.ampData = []
+        for extension in (0, 1, 2, 3):
+            calib.ampData.append(hdulist[extension].data)
+
+        calib.setMetadata(hdulist[0].header)
+
+        hdulist.close()
+        calib.log.info("Finished reading straylightData.")
+        return calib
 
     def evaluate(self, angle_start: Angle, angle_end: Optional[Angle] = None):
         """Get y-band background image array for a range of angles.
@@ -175,20 +192,19 @@ def evaluate(self, angle_start: Angle, angle_end: Optional[Angle] = None):
         ccd_img : `numpy.ndarray`
             Background data for this exposure.
         """
-        hdulist = fits.open(self._filename)
-        header = hdulist[0].header
+        header = self.getMetadata().toDict()
 
         # full-size ccd height & channel width
         ccd_h, ch_w = header["F_NAXIS2"], header["F_NAXIS1"]
         # saved data is compressed to 1/2**scale_level of the original size
         image_scale_level = header["WTLEVEL2"], header["WTLEVEL1"]
         angle_scale_level = header["WTLEVEL3"]
 
-        ccd_w = ch_w * len(hdulist)
+        ccd_w = ch_w * len(self.ampData)
         ccd_img = numpy.empty(shape=(ccd_h, ccd_w), dtype=numpy.float32)
 
-        for ch, hdu in enumerate(hdulist):
-            volume = _upscale_volume(hdu.data, angle_scale_level)
+        for ch, hdu in enumerate(self.ampData):
+            volume = _upscale_volume(hdu, angle_scale_level)
 
             if angle_end is None:
                 img = volume(angle_start.asDegrees())
@@ -207,18 +223,22 @@ def evaluate(self, angle_start: Angle, angle_end: Optional[Angle] = None):
 
 
 def _upscale_image(img, target_shape, level):
-    """
-    Upscale the given image to `target_shape` .
+    """Upscale the given image to `target_shape` .
 
-    @param img (numpy.array[][])
-        Compressed image. `img.shape` must agree
+    Parameters
+    ----------
+    img : `numpy.array`, (Nx, Ny)
+        Compressed image. ``img.shape`` must agree
         with waveletCompression.scaled_size(target_shape, level)
-    @param target_shape ((int, int))
+    target_shape : `tuple` [`int`, `int`]
         The shape of upscaled image, which is to be returned.
-    @param level (int or tuple of int)
+    level : `int` or `tuple` [`int`]
         Level of multiresolution analysis (or synthesis)
 
-    @return (numpy.array[][])
+    Returns
+    -------
+    resized : `numpy.array`, (Nu, Nv)
+        Upscaled image with the ``target_shape``.
     """
     h, w = waveletCompression.scaled_size(target_shape, level)
 
@@ -229,19 +249,23 @@ def _upscale_image(img, target_shape, level):
 
 
 def _upscale_volume(volume, level):
-    """
-    Upscale the given volume (= sequence of images) along the 0-th axis,
-    and return an instance of a interpolation object that interpolates
-    the 0-th axis. The 0-th axis is the instrument rotation.
+    """Upscale the given volume (= sequence of images) along the 0-th
+    axis, and return an instance of a interpolation object that
+    interpolates the 0-th axis. The 0-th axis is the instrument
+    rotation.
 
-    @param volume (numpy.array[][][])
+    Parameters
+    ----------
+    volume : `numpy.array`, (Nx, Ny, Nz)
         Sequence of images.
-    @param level (int)
+    level : `int`
         Level of multiresolution analysis along the 0-th axis.
 
-    @return (scipy.interpolate.CubicSpline)
-        You get a slice of the volume at a specific angle (in degrees)
-        by calling the returned value as `ret_value(angle)` .
+    interpolator : callable
+        An object that returns a slice of the volume at a specific
+        angle (in degrees), with one positional argument:
+
+        - ``angle``: The angle in degrees.
     """
     angles = 720
     _, h, w = volume.shape