From 7ae9d0cdffee8c0d85ec5487a62f16b32b469320 Mon Sep 17 00:00:00 2001 From: Ian Sullivan Date: Mon, 4 Nov 2024 20:57:53 -0800 Subject: [PATCH 1/3] Move common code to _prepareInputs method --- python/lsst/ip/diffim/detectAndMeasure.py | 42 +++++++++++++++++------ 1 file changed, 32 insertions(+), 10 deletions(-) diff --git a/python/lsst/ip/diffim/detectAndMeasure.py b/python/lsst/ip/diffim/detectAndMeasure.py index 881be05e..2f18ed86 100644 --- a/python/lsst/ip/diffim/detectAndMeasure.py +++ b/python/lsst/ip/diffim/detectAndMeasure.py @@ -189,6 +189,11 @@ class DetectAndMeasureConfig(pipeBase.PipelineTaskConfig, "base_PixelFlags_flag_saturatedCenterAll", ), ) + clearMaskPlanes = lsst.pex.config.ListField( + dtype=str, + doc="Mask planes to clear before running detection.", + default=("DETECTED", "DETECTED_NEGATIVE", "NOT_DEBLENDED", "STREAK"), + ) idGenerator = DetectorVisitIdGeneratorConfig.make_field() def setDefaults(self): @@ -339,13 +344,7 @@ def run(self, science, matchedTemplate, difference, if idFactory is None: idFactory = lsst.meas.base.IdGenerator().make_table_id_factory() - # Ensure that we start with an empty detection and deblended mask. - mask = difference.mask - clearMaskPlanes = ["DETECTED", "DETECTED_NEGATIVE", "NOT_DEBLENDED", "STREAK"] - for mp in clearMaskPlanes: - if mp not in mask.getMaskPlaneDict(): - mask.addMaskPlane(mp) - mask &= ~mask.getPlaneBitMask(clearMaskPlanes) + self._prepareInputs(difference) # Don't use the idFactory until after deblend+merge, so that we aren't # generating ids that just get thrown away (footprint merge doesn't @@ -365,6 +364,31 @@ def run(self, science, matchedTemplate, difference, positiveFootprints=positives, negativeFootprints=negatives) + def _prepareInputs(self, difference): + """Ensure that we start with an empty detection and deblended mask. + + Parameters + ---------- + difference : `lsst.afw.image.ExposureF` + The difference image that will be used for detecting diaSources. + The mask plane will be modified in place. + + Raises + ------ + lsst.pipe.base.UpstreamFailureNoWorkFound + If the PSF is not usable for measurement. + """ + # Check that we have a valid PSF now before we do more work + sigma = difference.psf.computeShape(difference.psf.getAveragePosition()).getDeterminantRadius() + if np.isnan(sigma): + raise pipeBase.UpstreamFailureNoWorkFound("Invalid PSF detected! PSF width evaluates to NaN.") + # Ensure that we start with an empty detection and deblended mask. + mask = difference.mask + for mp in self.config.clearMaskPlanes: + if mp not in mask.getMaskPlaneDict(): + mask.addMaskPlane(mp) + mask &= ~mask.getPlaneBitMask(self.config.clearMaskPlanes) + def processResults(self, science, matchedTemplate, difference, sources, idFactory, positiveFootprints=None, negativeFootprints=None,): """Measure and process the results of source detection. @@ -733,9 +757,7 @@ def run(self, science, matchedTemplate, difference, scoreExposure, if idFactory is None: idFactory = lsst.meas.base.IdGenerator().make_table_id_factory() - # Ensure that we start with an empty detection mask. - mask = scoreExposure.mask - mask &= ~(mask.getPlaneBitMask("DETECTED") | mask.getPlaneBitMask("DETECTED_NEGATIVE")) + self._prepareInputs(scoreExposure) # Don't use the idFactory until after deblend+merge, so that we aren't # generating ids that just get thrown away (footprint merge doesn't From 03400bccd527b3359ff6c91ecf17cbf461bc4dc4 Mon Sep 17 00:00:00 2001 From: Ian Sullivan Date: Tue, 5 Nov 2024 15:33:23 -0800 Subject: [PATCH 2/3] Do not decorrelate the diffim PSF for now. This is not implemented correctly, and raises occasional errors. --- python/lsst/ip/diffim/imageDecorrelation.py | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/python/lsst/ip/diffim/imageDecorrelation.py b/python/lsst/ip/diffim/imageDecorrelation.py index f85a4ae2..a2290482 100644 --- a/python/lsst/ip/diffim/imageDecorrelation.py +++ b/python/lsst/ip/diffim/imageDecorrelation.py @@ -274,7 +274,10 @@ def run(self, scienceExposure, templateExposure, subtractedExposure, psfMatching corr = self.computeDiffimCorrection(kArr, varianceMean, targetVarianceMean) correctedImage = self.computeCorrectedImage(corr.corrft, subtractedExposure.image.array) - correctedPsf = self.computeCorrectedDiffimPsf(corr.corrft, psfImg.array) + # TODO DM-47461: only decorrelate the PSF if it is calculated for the + # difference image. If it is taken directly from the science (or template) + # image the decorrelation calculation is incorrect. + # correctedPsf = self.computeCorrectedDiffimPsf(corr.corrft, psfImg.array) # The subtracted exposure variance plane is already correlated, we cannot propagate # it through another convolution; instead we need to use the uncorrelated originals @@ -300,7 +303,8 @@ def run(self, scienceExposure, templateExposure, subtractedExposure, psfMatching correctedVariance /= kSumSq subtractedExposure.image.array[...] = correctedImage # Allow for numpy type casting subtractedExposure.variance.array[...] = correctedVariance - subtractedExposure.setPsf(correctedPsf) + # TODO DM-47461 + # subtractedExposure.setPsf(correctedPsf) newVarMean = self.computeVarianceMean(subtractedExposure) self.log.info("Variance plane mean of corrected diffim: %.5e", newVarMean) From b991cbf3e96a3fc3c410942d4384aaad9da89c64 Mon Sep 17 00:00:00 2001 From: Ian Sullivan Date: Thu, 7 Nov 2024 17:26:36 -0800 Subject: [PATCH 3/3] Add unit test for bad PSF --- tests/test_detectAndMeasure.py | 33 ++++++++++++++++++++++++++++++++- 1 file changed, 32 insertions(+), 1 deletion(-) diff --git a/tests/test_detectAndMeasure.py b/tests/test_detectAndMeasure.py index da19000c..1b118162 100644 --- a/tests/test_detectAndMeasure.py +++ b/tests/test_detectAndMeasure.py @@ -22,9 +22,12 @@ import numpy as np import unittest +import lsst.afw.image as afwImage +import lsst.afw.math as afwMath import lsst.geom from lsst.ip.diffim import detectAndMeasure, subtractImages -from lsst.pipe.base import InvalidQuantumError +import lsst.meas.algorithms as measAlg +from lsst.pipe.base import InvalidQuantumError, UpstreamFailureNoWorkFound import lsst.utils.tests from utils import makeTestImage @@ -161,6 +164,34 @@ def test_detection_xy0(self): self.assertTrue(all(output.diaSources['id'] > 1000000000)) self.assertImagesEqual(subtractedMeasuredExposure.image, difference.image) + def test_raise_bad_psf(self): + """Detection should raise if the PSF width is NaN + """ + # Set up the simulated images + noiseLevel = 1. + staticSeed = 1 + fluxLevel = 500 + kwargs = {"seed": staticSeed, "psfSize": 2.4, "fluxLevel": fluxLevel, "x0": 12345, "y0": 67890} + science, sources = makeTestImage(noiseLevel=noiseLevel, noiseSeed=6, **kwargs) + matchedTemplate, _ = makeTestImage(noiseLevel=noiseLevel/4, noiseSeed=7, **kwargs) + difference = science.clone() + + psfShape = difference.getPsf().computeBBox(lsst.geom.Point2D(0, 0)).getDimensions() + psfcI = afwImage.ImageD(lsst.geom.Extent2I(psfShape[1], psfShape[0])) + psfNew = np.zeros(psfShape) + psfNew[0, :] = 1 + psfcI.array = psfNew + psfcK = afwMath.FixedKernel(psfcI) + psf = measAlg.KernelPsf(psfcK) + difference.setPsf(psf) + + # Configure the detection Task + detectionTask = self._setup_detection() + + # Run detection and check the results + with self.assertRaises(UpstreamFailureNoWorkFound): + detectionTask.run(science, matchedTemplate, difference) + def test_measurements_finite(self): """Measured fluxes and centroids should always be finite. """