DM-56872: Raise NoWorkFound for template errors with too few valid pixels

python/lsst/ip/diffim/getTemplate.py

@@ -53,7 +53,12 @@
 )
 from lsst.skymap import BaseSkyMap
 from lsst.ip.diffim.dcrModel import DcrModel
-from lsst.ip.diffim.utils import evaluateMeanPsfFwhm, getPsfFwhm, record_from_obs_info
+from lsst.ip.diffim.utils import (
+    checkTemplateIsSufficient,
+    evaluateMeanPsfFwhm,
+    getPsfFwhm,
+    record_from_obs_info,
+)
 from lsst.meas.algorithms import CoaddPsf, CoaddPsfConfig, SubtractBackgroundTask, ScaleVarianceTask
 from lsst.utils.timer import timeMethod
 
@@ -340,6 +345,12 @@ class GetTemplateConfig(
         dtype=float,
         default=0.05,
     )
+    minTemplateFractionForExpectedSuccess = pexConfig.Field(
+        doc="Raise NoWorkFound if the template PSF cannot be evaluated and less than this fraction"
+            " of the template pixels have data. With more coverage, the PSF error is raised.",
+        dtype=float,
+        default=0.1,
+    )
     varianceBackground = pexConfig.ConfigurableField(
         target=SubtractBackgroundTask,
         doc="Task to estimate the background variance.",
@@ -732,7 +743,12 @@ def run(self, *, coaddExposureHandles, bbox, wcs, dataIds, physical_filter, visi
         Raises
         ------
         NoWorkFound
-            If no coadds are found with sufficient un-masked pixels.
+            If no coadds are found with sufficient un-masked pixels, or if the
+            template PSF cannot be evaluated and too few of the template pixels
+            have data.
+        lsst.afw.detection.InvalidPsfError
+            If the template PSF cannot be evaluated and the template has
+            enough coverage that it was expected to succeed.
         RuntimeError
             If the template is not known to be from the science image's
             instrument and no coadd records the plate scale of its input
@@ -841,7 +857,19 @@ def run(self, *, coaddExposureHandles, bbox, wcs, dataIds, physical_filter, visi
             photoCalib = identity
         template.setPhotoCalib(photoCalib)
         template.setPsf(self._makePsf(template, catalog, wcs))
-        template.setPsf(self._makeGaussianPsf(template))
+        try:
+            template.setPsf(self._makeGaussianPsf(template))
+        except afwDetection.InvalidPsfError as e:
+            # Raise NoWorkFound if the template coverage is insufficient.
+            checkTemplateIsSufficient(
+                template, self.log,
+                requiredTemplateFraction=self.config.minTemplateFractionForExpectedSuccess,
+                exceptionMessage="Template coverage lower than expected to succeed."
+                f" Failure is tolerable: {e}"
+            )
+            # checkTemplateIsSufficient did not raise NoWorkFound, so raise
+            # the original exception.
+            raise
 
         # Record the input coadd patches as the template's coadd inputs.
         coaddInputs = afwImage.CoaddInputs(afwTable.ExposureTable.makeMinimalSchema(), self.schema)
@@ -1304,7 +1332,7 @@ def _makeGaussianPsf(self, template):
 
         Raises
         ------
-        ValueError
+        lsst.afw.detection.InvalidPsfError
             Raised if the PSF cannot be evaluated at any grid position.
         """
         # The CoaddPsf is not defined where the template has no inputs, which

python/lsst/ip/diffim/subtractImages.py

@@ -34,7 +34,7 @@
 from lsst.ip.diffim.utils import (evaluateMeanPsfFwhm, getPsfFwhm,
                                   computeDifferenceImageMetrics, computePSFNoiseEquivalentArea,
                                   checkMask, record_from_obs_info, setSourceFootprints,
-                                  get_difference_image_provenance)
+                                  get_difference_image_provenance, checkTemplateIsSufficient)
 from lsst.meas.algorithms import ScaleVarianceTask, ScienceSourceSelectorTask
 import lsst.pex.config
 import lsst.pipe.base
@@ -749,10 +749,12 @@ def runMakeKernel(self, template, science, sources=None, convolveTemplate=True,
         except (RuntimeError, lsst.pex.exceptions.Exception) as e:
             self.log.warning("Failed to match template. Checking coverage")
             #  Raise NoWorkFound if template fraction is insufficient
-            checkTemplateIsSufficient(template[science.getBBox()], science, self.log,
-                                      self.config.minTemplateFractionForExpectedSuccess,
-                                      exceptionMessage="Template coverage lower than expected to succeed."
-                                      f" Failure is tolerable: {e}")
+            checkTemplateIsSufficient(
+                template[science.getBBox()], self.log, scienceExposure=science,
+                requiredTemplateFraction=self.config.minTemplateFractionForExpectedSuccess,
+                exceptionMessage="Template coverage lower than expected to succeed."
+                f" Failure is tolerable: {e}"
+            )
             #  checkTemplateIsSufficient did not raise NoWorkFound, so raise original exception
             raise e
 
@@ -1192,7 +1194,7 @@ def _prepareInputs(self, template, science, visitSummary=None):
         if visitSummary is not None:
             self._applyExternalCalibrations(science, visitSummary=visitSummary)
         templateCoverageFraction = checkTemplateIsSufficient(
-            template[science.getBBox()], science, self.log,
+            template[science.getBBox()], self.log, scienceExposure=science,
             requiredTemplateFraction=self.config.requiredTemplateFraction,
             exceptionMessage="Not attempting subtraction. To force subtraction,"
             " set config requiredTemplateFraction=0"
@@ -1432,10 +1434,12 @@ def run(self, template, science, sources, visitSummary=None):
         except (RuntimeError, lsst.pex.exceptions.Exception) as e:
             self.log.warning("Failed to match template. Checking coverage")
             #  Raise NoWorkFound if template fraction is insufficient
-            checkTemplateIsSufficient(template[science.getBBox()], science, self.log,
-                                      self.config.minTemplateFractionForExpectedSuccess,
-                                      exceptionMessage="Template coverage lower than expected to succeed."
-                                      f" Failure is tolerable: {e}")
+            checkTemplateIsSufficient(
+                template[science.getBBox()], self.log, scienceExposure=science,
+                requiredTemplateFraction=self.config.minTemplateFractionForExpectedSuccess,
+                exceptionMessage="Template coverage lower than expected to succeed."
+                f" Failure is tolerable: {e}"
+            )
             #  checkTemplateIsSufficient did not raise NoWorkFound, so raise original exception
             raise e
 
@@ -1583,53 +1587,6 @@ def runPreconvolve(self, template, science, matchedScience, kernelSources, preCo
                                      kernelSources=kernelSources)
 
 
-def checkTemplateIsSufficient(templateExposure, scienceExposure, logger, requiredTemplateFraction=0.,
-                              exceptionMessage=""):
-    """Raise NoWorkFound if template coverage < requiredTemplateFraction
-
-    Parameters
-    ----------
-    templateExposure : `lsst.afw.image.ExposureF`
-        The template exposure to check
-    logger : `logging.Logger`
-        Logger for printing output.
-    requiredTemplateFraction : `float`, optional
-        Fraction of pixels of the science image required to have coverage
-        in the template.
-    exceptionMessage : `str`, optional
-        Message to include in the exception raised if the template coverage
-        is insufficient.
-
-    Returns
-    -------
-    templateCoverageFraction: `float`
-        Fraction of pixels in the template with data.
-
-    Raises
-    ------
-    lsst.pipe.base.NoWorkFound
-        Raised if fraction of good pixels, defined as not having NO_DATA
-        set, is less than the requiredTemplateFraction
-    """
-    # Count the number of pixels with the NO_DATA mask bit set
-    # counting NaN pixels is insufficient because pixels without data are often intepolated over)
-    noTemplate = templateExposure.mask.array & templateExposure.mask.getPlaneBitMask('NO_DATA')
-    # Also need to account for missing data in the science image,
-    # because template coverage there doesn't help
-    noScience = scienceExposure.mask.array & scienceExposure.mask.getPlaneBitMask('NO_DATA')
-    pixNoData = np.count_nonzero(noTemplate | noScience)
-    pixGood = templateExposure.getBBox().getArea() - pixNoData
-    templateCoverageFraction = pixGood/templateExposure.getBBox().getArea()
-    logger.info("template has %d good pixels (%.1f%%)", pixGood, 100*templateCoverageFraction)
-
-    if templateCoverageFraction < requiredTemplateFraction:
-        message = ("Insufficient Template Coverage. (%.1f%% < %.1f%%)" % (
-                   100*templateCoverageFraction,
-                   100*requiredTemplateFraction))
-        raise lsst.pipe.base.NoWorkFound(message + " " + exceptionMessage)
-    return templateCoverageFraction
-
-
 def _subtractImages(science, template, backgroundModel=None):
     """Subtract template from science, propagating relevant metadata.
 

python/lsst/ip/diffim/utils.py

@@ -25,7 +25,7 @@
 __all__ = ["evaluateMeanPsfFwhm", "getPsfFwhm", "getKernelCenterDisplacement",
            "computeDifferenceImageMetrics", "checkMask", "setSourceFootprints",
            "record_from_obs_info", "get_difference_image_provenance",
-           ]
+           "checkTemplateIsSufficient"]
 
 import itertools
 import os
@@ -266,7 +266,7 @@ def evaluateMeanPsfFwhm(exposure: afwImage.Exposure,
 
     Raises
     ------
-    ValueError
+    lsst.afw.detection.InvalidPsfError
         Raised if the PSF cannot be computed at any of the grid points.
 
     See Also
@@ -298,7 +298,7 @@ def evaluateMeanPsfFwhm(exposure: afwImage.Exposure,
         width.append(fwhm)
 
     if not width:
-        raise ValueError("Unable to compute PSF FWHM at any position on the exposure.")
+        raise afwDetection.InvalidPsfError("Unable to compute PSF FWHM at any position on the exposure.")
 
     return np.nanmean(width)
 
@@ -633,3 +633,55 @@ def setSourceFootprints(sources, kernelSize):
         boxFootprint.addPeak(peak.getFx(), peak.getFy(), peak.getPeakValue())
         source.setFootprint(boxFootprint)
     return sources
+
+
+def checkTemplateIsSufficient(templateExposure, logger, scienceExposure=None, requiredTemplateFraction=0.,
+                              exceptionMessage=""):
+    """Raise NoWorkFound if template coverage < requiredTemplateFraction
+
+    Parameters
+    ----------
+    templateExposure : `lsst.afw.image.ExposureF`
+        The template exposure to check
+    logger : `logging.Logger`
+        Logger for printing output.
+    scienceExposure : `lsst.afw.image.ExposureF`, optional
+        The science image, used to constrain the template pixels.
+    requiredTemplateFraction : `float`, optional
+        Fraction of pixels of the science image required to have coverage
+        in the template.
+    exceptionMessage : `str`, optional
+        Message to include in the exception raised if the template coverage
+        is insufficient.
+
+    Returns
+    -------
+    templateCoverageFraction: `float`
+        Fraction of pixels in the template with data.
+
+    Raises
+    ------
+    lsst.pipe.base.NoWorkFound
+        Raised if fraction of good pixels, defined as not having NO_DATA
+        set, is less than the requiredTemplateFraction
+    """
+    # Count the pixels with the NO_DATA mask bit set. Counting NaN pixels is
+    # insufficient because pixels without data are often interpolated over.
+    noTemplate = templateExposure.mask.array & templateExposure.mask.getPlaneBitMask('NO_DATA')
+    # Also need to account for missing data in the science image,
+    # because template coverage there doesn't help
+    if scienceExposure is not None:
+        noScience = scienceExposure.mask.array & scienceExposure.mask.getPlaneBitMask('NO_DATA')
+    else:
+        noScience = False
+    pixNoData = np.count_nonzero(noTemplate | noScience)
+    pixGood = templateExposure.getBBox().getArea() - pixNoData
+    templateCoverageFraction = pixGood/templateExposure.getBBox().getArea()
+    logger.info("template has %d good pixels (%.1f%%)", pixGood, 100*templateCoverageFraction)
+
+    if templateCoverageFraction < requiredTemplateFraction:
+        message = ("Insufficient Template Coverage. (%.1f%% < %.1f%%)" % (
+                   100*templateCoverageFraction,
+                   100*requiredTemplateFraction))
+        raise lsst.pipe.base.NoWorkFound(message + " " + exceptionMessage)
+    return templateCoverageFraction

tests/test_getTemplate.py

@@ -403,6 +403,30 @@ def testGaussianPsfGridFallback(self):
         self.assertFloatsAlmostEqual(psf.getSigma(), expected)
         self.assertFloatsAlmostEqual(task.metadata["templatePsfSigma"], expected)
 
+    @lsst.utils.tests.methodParameters(
+        minFraction=[0.1, 0.9],
+        expected=[lsst.afw.detection.InvalidPsfError, pipeBase.NoWorkFound],
+    )
+    def testGaussianPsfFailure(self, minFraction=None, expected=None):
+        """A template PSF that cannot be evaluated raises NoWorkFound if the
+        template coverage is below the configured minimum, and the PSF error
+        otherwise.
+        """
+        # Drop a patch so that only about half of the template has data.
+        self.patches[0].pop(1)
+        box = lsst.geom.Box2I(lsst.geom.Point2I(0, 0), lsst.geom.Point2I(180, 180))
+        config = self._makeConfig()
+        config.minTemplateFractionForExpectedSuccess = minFraction
+        task = lsst.ip.diffim.GetTemplateTask(config=config)
+        error = lsst.afw.detection.InvalidPsfError("Unable to compute PSF FWHM.")
+        with unittest.mock.patch.object(task, "_makeGaussianPsf", side_effect=error):
+            with self.assertRaisesRegex(expected, "Unable to compute PSF FWHM"):
+                task.run(coaddExposureHandles=self.patches,
+                         bbox=box,
+                         wcs=self.exposure.wcs,
+                         dataIds=self.dataIds,
+                         physical_filter="a_test")
+
     def testRunNoTemplate(self):
         """A bounding box that doesn't overlap the patches will raise.
         """

tests/test_utils.py

@@ -20,14 +20,23 @@
 # along with this program.  If not, see <https://www.gnu.org/licenses/>.
 
 
+import logging
 import requests
 from unittest import mock
 
 import lsst.afw.cameraGeom.testUtils
+import lsst.afw.detection
 import lsst.daf.butler
 from lsst.images import obs_info_from_legacy
 from lsst.images.tests import get_dp2_exposure_record
-from lsst.ip.diffim.utils import populate_sattle_visit_cache, record_from_obs_info
+from lsst.ip.diffim.utils import (
+    checkTemplateIsSufficient,
+    evaluateMeanPsfFwhm,
+    populate_sattle_visit_cache,
+    record_from_obs_info,
+)
+from lsst.pex.exceptions import InvalidParameterError
+import lsst.pipe.base
 import lsst.utils.tests
 
 from test_detectAndMeasure import makeVisitInfo, MockResponse
@@ -100,3 +109,29 @@ def test_populate_sattle_raises(self):
     def test_raise_on_even_kernel(self):
         with self.assertRaises(ValueError):
             makeTestImage(kernelSize=32)
+
+    def test_evaluateMeanPsfFwhm_raises(self):
+        """A PSF that cannot be evaluated at any grid point is invalid."""
+        exposure, _ = makeTestImage(nSrc=1, xSize=64, ySize=64)
+        error = InvalidParameterError("No inputs at this position.")
+        with mock.patch("lsst.ip.diffim.utils.getPsfFwhm", side_effect=error):
+            with self.assertRaisesRegex(lsst.afw.detection.InvalidPsfError, "at any position"):
+                evaluateMeanPsfFwhm(exposure, fwhmExposureBuffer=0.05, fwhmExposureGrid=3)
+
+    def test_checkTemplateIsSufficient(self):
+        """The coverage is the fraction of pixels with data in the template
+        and, if supplied, also in the science image.
+        """
+        template, _ = makeTestImage(nSrc=1, xSize=64, ySize=64, clearEdgeMask=True)
+        science = template.clone()
+        noData = template.mask.getPlaneBitMask("NO_DATA")
+        template.mask.array[:32, :] |= noData
+        science.mask.array[:, :32] |= noData
+        logger = logging.getLogger(__name__)
+
+        self.assertEqual(checkTemplateIsSufficient(template, logger), 0.5)
+        self.assertEqual(checkTemplateIsSufficient(template, logger, scienceExposure=science), 0.25)
+        self.assertEqual(checkTemplateIsSufficient(template, logger, requiredTemplateFraction=0.5), 0.5)
+        with self.assertRaisesRegex(lsst.pipe.base.NoWorkFound, "Insufficient Template Coverage"):
+            checkTemplateIsSufficient(template, logger, scienceExposure=science,
+                                      requiredTemplateFraction=0.5)