From b6da15f71bcf562fb5f425a178fd3a6830c78dc0 Mon Sep 17 00:00:00 2001 From: Arun Kannawadi Date: Thu, 20 Feb 2025 12:47:49 -0500 Subject: [PATCH 1/6] Move up clearing NaNs from weight Having NaN values in the weight that multiplies images is bad, since a single bad image can ruin a pixel regardless of other inputs. Since we are setting its corresponding weight to be zero, we kill its image (and variance) values also by multiplying by zero. --- python/lsst/ip/diffim/getTemplate.py | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/python/lsst/ip/diffim/getTemplate.py b/python/lsst/ip/diffim/getTemplate.py index 04e9d5ce6..379ab1bfd 100644 --- a/python/lsst/ip/diffim/getTemplate.py +++ b/python/lsst/ip/diffim/getTemplate.py @@ -406,16 +406,17 @@ def _merge(self, maskedImages, bbox, wcs): maskedImage.variance.array[bad] = 0.0 # Reset mask, too, since these pixels don't contribute to sum. maskedImage.mask.array[bad] = 0 + # Clear the NaNs to ensure that areas missing from this input are + # masked with NO_DATA after the loop. + weight.array[np.isnan(weight.array)] = 0 # Cannot use `merged.maskedImage *= weight` because that operator # multiplies the variance by the weight twice; in this case # `weight` are the exact values we want to scale by. maskedImage.image *= weight maskedImage.variance *= weight merged.maskedImage[maskedImage.getBBox()] += maskedImage - # Clear the NaNs to ensure that areas missing from this input are - # masked with NO_DATA after the loop. - weight.array[np.isnan(weight.array)] = 0 weights[maskedImage.getBBox()] += weight + # Cannot use `merged.maskedImage /= weights` because that operator # divides the variance by the weight twice; in this case `weights` are # the exact values we want to scale by. From 428322f93fb04f2a6c1aaa838d1f5cd75d5af5ea Mon Sep 17 00:00:00 2001 From: Arun Kannawadi Date: Thu, 20 Feb 2025 12:49:33 -0500 Subject: [PATCH 2/6] Weed out infinities In case the variance plane has zeros resulting in inf values in weight, that needs to be cleared out as well. Instead of mapping only NaN values to zero, we mapping non-finite values to zero more generally. --- python/lsst/ip/diffim/getTemplate.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/python/lsst/ip/diffim/getTemplate.py b/python/lsst/ip/diffim/getTemplate.py index 379ab1bfd..ede2ea50b 100644 --- a/python/lsst/ip/diffim/getTemplate.py +++ b/python/lsst/ip/diffim/getTemplate.py @@ -407,8 +407,9 @@ def _merge(self, maskedImages, bbox, wcs): # Reset mask, too, since these pixels don't contribute to sum. maskedImage.mask.array[bad] = 0 # Clear the NaNs to ensure that areas missing from this input are - # masked with NO_DATA after the loop. - weight.array[np.isnan(weight.array)] = 0 + # masked with NO_DATA after the loop. Limiting to finite values + # can also handle the case where the input variance was zeros. + weight.array[~np.isfinite(weight.array)] = 0 # Cannot use `merged.maskedImage *= weight` because that operator # multiplies the variance by the weight twice; in this case # `weight` are the exact values we want to scale by. From e60d4feb8958d5b453b06e11c93e734e6d087ed6 Mon Sep 17 00:00:00 2001 From: Arun Kannawadi Date: Thu, 20 Feb 2025 12:50:53 -0500 Subject: [PATCH 3/6] Guard against division by zero It is inconsistent to set quantities to zero, and not guard against division by the same quantities. Right before the division, we map zeros to inf, so the image and variance values of those pixels are zeros. We also set NO_DATA mask bit when that happens. --- python/lsst/ip/diffim/getTemplate.py | 4 +++- tests/test_getTemplate.py | 6 +++--- 2 files changed, 6 insertions(+), 4 deletions(-) diff --git a/python/lsst/ip/diffim/getTemplate.py b/python/lsst/ip/diffim/getTemplate.py index ede2ea50b..42491a4de 100644 --- a/python/lsst/ip/diffim/getTemplate.py +++ b/python/lsst/ip/diffim/getTemplate.py @@ -418,12 +418,14 @@ def _merge(self, maskedImages, bbox, wcs): merged.maskedImage[maskedImage.getBBox()] += maskedImage weights[maskedImage.getBBox()] += weight + bad = np.isnan(weights.array) | (weights.array == 0) + weights.array[bad] = np.inf # Cannot use `merged.maskedImage /= weights` because that operator # divides the variance by the weight twice; in this case `weights` are # the exact values we want to scale by. merged.image /= weights merged.variance /= weights - merged.mask.array |= merged.mask.getPlaneBitMask("NO_DATA") * (weights.array == 0) + merged.mask.array |= merged.mask.getPlaneBitMask("NO_DATA") * bad return merged diff --git a/tests/test_getTemplate.py b/tests/test_getTemplate.py index a9b00e9b1..b8c728a76 100644 --- a/tests/test_getTemplate.py +++ b/tests/test_getTemplate.py @@ -282,9 +282,9 @@ def testMissingPatches(self): dataIds=self.dataIds, physical_filter="a_test") no_data = (result.template.mask.array & result.template.mask.getPlaneBitMask("NO_DATA")) != 0 - self.assertTrue(all(np.isnan(result.template.image.array[no_data]))) - self.assertTrue(all(np.isnan(result.template.variance.array[no_data]))) - self.assertEqual(no_data.sum(), 21548) + self.assertTrue(np.isfinite(result.template.image.array).all()) + self.assertTrue(np.isfinite(result.template.variance.array).all()) + self.assertEqual(no_data.sum(), 20990) def setup_module(module): From c219978ef893c31744eefc757a56604fa95d35b2 Mon Sep 17 00:00:00 2001 From: Arun Kannawadi Date: Thu, 20 Feb 2025 12:53:21 -0500 Subject: [PATCH 4/6] Add a unit test for inputs with NaN This unit test mimics a condition that triggered DM-48920. --- tests/test_getTemplate.py | 28 ++++++++++++++++++++++++++++ 1 file changed, 28 insertions(+) diff --git a/tests/test_getTemplate.py b/tests/test_getTemplate.py index b8c728a76..1668c8b13 100644 --- a/tests/test_getTemplate.py +++ b/tests/test_getTemplate.py @@ -286,6 +286,34 @@ def testMissingPatches(self): self.assertTrue(np.isfinite(result.template.variance.array).all()) self.assertEqual(no_data.sum(), 20990) + @lsst.utils.tests.methodParameters( + box=[ + lsst.geom.Box2I(lsst.geom.Point2I(0, 0), lsst.geom.Point2I(180, 180)), + lsst.geom.Box2I(lsst.geom.Point2I(200, 200), lsst.geom.Point2I(600, 600)), + ], + nInput=[8, 16], + ) + def testNanInputs(self, box=None, nInput=None): + """Test that the template has finite values when some of the input pixels have NaN as variance.""" + for tract, patchCoadds in self.patches.items(): + for patchCoadd in patchCoadds: + bbox = lsst.geom.Box2I() + bbox.include(lsst.geom.Point2I(patchCoadd.getBBox().getCenter())) + bbox.grow(3) + patchCoadd.variance[bbox].array *= np.nan + + box = lsst.geom.Box2I(lsst.geom.Point2I(200, 200), lsst.geom.Point2I(600, 600)) + task = lsst.ip.diffim.GetTemplateTask() + result = task.run(coaddExposures=self.patches, + bbox=lsst.geom.Box2I(box), + wcs=self.exposure.wcs, + dataIds=self.dataIds, + physical_filter="a_test") + self._checkMetadata(result.template, task.config, box, self.exposure.wcs, 16) + # We just check that the pixel values are all finite. We cannot check that pixel values + # in the template are closer to the original anymore. + self.assertTrue(np.isfinite(result.template.image.array).all()) + def setup_module(module): lsst.utils.tests.init() From 89c90a54b6708ce5fb2320ea2e958f5371f4f4af Mon Sep 17 00:00:00 2001 From: Arun Kannawadi Date: Thu, 20 Feb 2025 15:45:19 -0500 Subject: [PATCH 5/6] Declare _merge as a staticmethod --- python/lsst/ip/diffim/getTemplate.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/python/lsst/ip/diffim/getTemplate.py b/python/lsst/ip/diffim/getTemplate.py index 42491a4de..82873b6b8 100644 --- a/python/lsst/ip/diffim/getTemplate.py +++ b/python/lsst/ip/diffim/getTemplate.py @@ -375,7 +375,8 @@ def _makeExposureCatalog(self, exposures, dataIds): return images, catalog, totalBox - def _merge(self, maskedImages, bbox, wcs): + @staticmethod + def _merge(maskedImages, bbox, wcs): """Merge the images that came from one tract into one larger image, ignoring NaN pixels and non-finite variance pixels from individual exposures. From 40dcd7f7159c4264acf4e67e71c160ffe43801f2 Mon Sep 17 00:00:00 2001 From: Ian Sullivan Date: Mon, 3 Mar 2025 13:43:52 -0800 Subject: [PATCH 6/6] Reconcile changes with DM-49268 These additional changes address more of the divide by zero and NaN warnings --- python/lsst/ip/diffim/getTemplate.py | 29 ++++++++++++++-------------- 1 file changed, 15 insertions(+), 14 deletions(-) diff --git a/python/lsst/ip/diffim/getTemplate.py b/python/lsst/ip/diffim/getTemplate.py index 82873b6b8..1d1d12523 100644 --- a/python/lsst/ip/diffim/getTemplate.py +++ b/python/lsst/ip/diffim/getTemplate.py @@ -399,18 +399,17 @@ def _merge(maskedImages, bbox, wcs): merged = afwImage.ExposureF(bbox, wcs) weights = afwImage.ImageF(bbox) for maskedImage in maskedImages: - weight = afwImage.ImageF(maskedImage.variance.array**(-0.5)) - bad = np.isnan(maskedImage.image.array) | ~np.isfinite(maskedImage.variance.array) + # Catch both zero-value and NaN variance plane pixels + good = maskedImage.variance.array > 0 + weight = afwImage.ImageF(maskedImage.getBBox()) + weight.array[good] = maskedImage.variance.array[good]**(-0.5) + bad = np.isnan(maskedImage.image.array) | ~good # Note that modifying the patch MaskedImage in place is fine; # we're throwing it away at the end anyway. maskedImage.image.array[bad] = 0.0 maskedImage.variance.array[bad] = 0.0 # Reset mask, too, since these pixels don't contribute to sum. maskedImage.mask.array[bad] = 0 - # Clear the NaNs to ensure that areas missing from this input are - # masked with NO_DATA after the loop. Limiting to finite values - # can also handle the case where the input variance was zeros. - weight.array[~np.isfinite(weight.array)] = 0 # Cannot use `merged.maskedImage *= weight` because that operator # multiplies the variance by the weight twice; in this case # `weight` are the exact values we want to scale by. @@ -419,14 +418,16 @@ def _merge(maskedImages, bbox, wcs): merged.maskedImage[maskedImage.getBBox()] += maskedImage weights[maskedImage.getBBox()] += weight - bad = np.isnan(weights.array) | (weights.array == 0) - weights.array[bad] = np.inf - # Cannot use `merged.maskedImage /= weights` because that operator - # divides the variance by the weight twice; in this case `weights` are - # the exact values we want to scale by. - merged.image /= weights - merged.variance /= weights - merged.mask.array |= merged.mask.getPlaneBitMask("NO_DATA") * bad + inverseWeights = np.zeros_like(weights.array) + good = weights.array > 0 + inverseWeights[good] = 1/weights.array[good] + + # Cannot use `merged.maskedImage *= inverseWeights` because that + # operator divides the variance by the weight twice; in this case + # `inverseWeights` are the exact values we want to scale by. + merged.image.array *= inverseWeights + merged.variance.array *= inverseWeights + merged.mask.array |= merged.mask.getPlaneBitMask("NO_DATA") * (inverseWeights == 0) return merged