From 7b399bd07e062e37c431eebad8e5d5d6a89c250b Mon Sep 17 00:00:00 2001 From: Merlin Fisher-Levine Date: Thu, 29 Sep 2022 13:13:43 -0700 Subject: [PATCH 1/7] Add utility function for summing footprint fluxes --- python/lsst/summit/utils/utils.py | 76 +++++++++++++++++++++++++++++++ 1 file changed, 76 insertions(+) diff --git a/python/lsst/summit/utils/utils.py b/python/lsst/summit/utils/utils.py index 150bb8e8..217d5924 100644 --- a/python/lsst/summit/utils/utils.py +++ b/python/lsst/summit/utils/utils.py @@ -20,11 +20,13 @@ # along with this program. If not, see . import os +from typing import Iterable import numpy as np import logging from scipy.ndimage import gaussian_filter import lsst.afw.image as afwImage import lsst.afw.detection as afwDetect +from lsst.afw.detection import Footprint, FootprintSet import lsst.afw.math as afwMath import lsst.daf.base as dafBase import lsst.geom as geom @@ -308,6 +310,80 @@ def detectObjectsInExp(exp, nSigma=10, nPixMin=10, grow=0): return footPrintSet +def fluxesFromFootprints(footprints, parentImage, subtractImageMedian=False): + """Calculate the flux from a set of footprints, given the parent image, + optionally subtracting the whole-image median from each pixel as a very + rough background subtraction. + + Parameters + ---------- + footprints : `lsst.afw.detection.FootprintSet` or + `lsst.afw.detection.Footprint` or + `iterable` of `lsst.afw.detection.Footprint` + The footprints to measure. + parentImage : `lsst.afw.image.Exposure` + The parent exposure. + subtractImageMedian : `bool`, optional + Subtract a whole-image median from each pixel in the footprint when + summing as a very crude background subtraction. Does not change the + original image. + + Returns + ------- + fluxes : `list` of `float` + The fluxes for each footprint. + + Raises + ------ + ValueError : raise for unsupported types. + """ + median = 0 + if subtractImageMedian: + median = np.nanmedian(parentImage.image.array) + + # poor person's single dispatch + badTypeMsg = ("This function works with FootprintSets, single Footprints, and iterables of Footprints. " + f"Got {type(footprints)}: {footprints}") + if isinstance(footprints, FootprintSet): + fps = footprints.getFootprints() + elif isinstance(footprints, Iterable): + if not isinstance(footprints[0], Footprint): + raise ValueError(badTypeMsg) + fps = footprints + elif isinstance(footprints, Footprint): + fps = [footprints] + else: + raise ValueError(badTypeMsg) + + return [fluxFromFootprint(fp, parentImage, backgroundValue=median) for fp in fps] + + +def fluxFromFootprint(footprint, parentImage, backgroundValue=0): + """Calculate the flux from a footprint, given the parent image, optionally + subtracting a single value from each pixel as a very rough background + subtraction, e.g. the image median. + + Parameters + ---------- + footprint : `lsst.afw.detection.Footprint` + The footprint to measure. + parentImage : `lsst.afw.image.Exposure` + The parent exposure. + backgroundValue : `bool`, optional + The value to subtract from each pixel in the footprint when summing + as a very crude background subtraction. Does not change the original + image. + + Returns + ------- + flux : `float` + The flux in the footprint + """ + if backgroundValue: # only do the subtraction if non-zero for speed + return footprint.getSpans().flatten(parentImage.image.array - backgroundValue).sum() + return footprint.getSpans().flatten(parentImage.image.array).sum() + + def humanNameForCelestialObject(objName): """Returns a list of all human names for obj, or [] if none are found. From 1a84497b04caf83008fd04d3f725741fe2c5023e Mon Sep 17 00:00:00 2001 From: Merlin Fisher-Levine Date: Mon, 3 Oct 2022 18:47:38 -0700 Subject: [PATCH 2/7] Add tests for flux calculating functions --- tests/test_utils.py | 47 +++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 47 insertions(+) diff --git a/tests/test_utils.py b/tests/test_utils.py index 8d63cc82..98e7c737 100644 --- a/tests/test_utils.py +++ b/tests/test_utils.py @@ -23,6 +23,7 @@ import copy import itertools +from typing import Iterable import unittest import astropy.time @@ -35,14 +36,21 @@ from lsst.obs.base import createInitialSkyWcsFromBoresight from lsst.obs.base.makeRawVisitInfoViaObsInfo import MakeRawVisitInfoViaObsInfo from lsst.obs.lsst import Latiss + from lsst.summit.utils.utils import (getExpPositionOffset, getFieldNameAndTileNumber, getAirmassSeeingCorrection, getFilterSeeingCorrection, quickSmooth, getQuantiles, + fluxesFromFootprints, + detectObjectsInExp, ) +from lsst.summit.utils.bestEffort import BestEffortIsr +from lsst.summit.utils.butlerUtils import makeDefaultLatissButler + from lsst.obs.lsst.translators.latiss import AUXTEL_LOCATION +from lsst.afw.detection import FootprintSet class ExpSkyPositionOffsetTestCase(lsst.utils.tests.TestCase): @@ -207,6 +215,45 @@ def test_quantiles(self): np.testing.assert_almost_equal(edges1, edges2, decimal=decimal) +class ImageBasedTestCase(lsst.utils.tests.TestCase): + """A test case for testing sky position offsets for exposures.""" + + @classmethod + def setUpClass(cls): + try: + cls.butler = makeDefaultLatissButler() + except FileNotFoundError: + raise unittest.SkipTest("Skipping tests that require the LATISS butler repo.") + + cls.dataId = {'day_obs': 20200315, 'seq_num': 120, 'detector': 0} + cls.bestEffort = BestEffortIsr() # should always succeed if makeDefaultLatissButler works + cls.exp = cls.bestEffort.getExposure(cls.dataId) + + def test_fluxFromFootprint(self): + footPrintSet = detectObjectsInExp(self.exp) + allFootprints = footPrintSet.getFootprints() + self.assertGreaterEqual(len(allFootprints), 3) + singleFootprint = allFootprints[0] + twoFootprints = allFootprints[0:2] + + # check it can accept a footprintSet, and single and iterables of + # footprints + for obj in (footPrintSet, singleFootprint, twoFootprints): + fluxes = fluxesFromFootprints(obj, self.exp) + if isinstance(obj, FootprintSet): + expectedLength = len(footPrintSet.getFootprints()) + self.assertEqual(len(fluxes), expectedLength) # always one flux per footprint + + self.assertIsInstance(fluxes, Iterable) + self.assertIsInstance(fluxes[0], np.floating) # always returns floats + + # ensure that subtracting the image median from fluxes leave image + # functionally unchanged + oldExpArray = copy.deepcopy(self.exp.image.array) + fluxes = fluxesFromFootprints(obj, self.exp, subtractImageMedian=True) + self.assertTrue(np.allclose(self.exp.image.array, oldExpArray)) + + class TestMemory(lsst.utils.tests.MemoryTestCase): pass From 6771e1568e578cd0e4e44df39866e76d2bafb9f0 Mon Sep 17 00:00:00 2001 From: Merlin Fisher-Levine Date: Tue, 25 Oct 2022 07:35:19 -0700 Subject: [PATCH 3/7] Change raise from ValueError to TypeError --- python/lsst/summit/utils/utils.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/python/lsst/summit/utils/utils.py b/python/lsst/summit/utils/utils.py index 217d5924..53ffe80e 100644 --- a/python/lsst/summit/utils/utils.py +++ b/python/lsst/summit/utils/utils.py @@ -335,7 +335,7 @@ def fluxesFromFootprints(footprints, parentImage, subtractImageMedian=False): Raises ------ - ValueError : raise for unsupported types. + TypeError : raise for unsupported types. """ median = 0 if subtractImageMedian: @@ -348,12 +348,12 @@ def fluxesFromFootprints(footprints, parentImage, subtractImageMedian=False): fps = footprints.getFootprints() elif isinstance(footprints, Iterable): if not isinstance(footprints[0], Footprint): - raise ValueError(badTypeMsg) + raise TypeError(badTypeMsg) fps = footprints elif isinstance(footprints, Footprint): fps = [footprints] else: - raise ValueError(badTypeMsg) + raise TypeError(badTypeMsg) return [fluxFromFootprint(fp, parentImage, backgroundValue=median) for fp in fps] From 0820479ff12c2028998d15632a36df931c2493cb Mon Sep 17 00:00:00 2001 From: Merlin Fisher-Levine Date: Tue, 25 Oct 2022 07:43:16 -0700 Subject: [PATCH 4/7] Remove use of fps in single dispatch --- python/lsst/summit/utils/utils.py | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/python/lsst/summit/utils/utils.py b/python/lsst/summit/utils/utils.py index 53ffe80e..00440de4 100644 --- a/python/lsst/summit/utils/utils.py +++ b/python/lsst/summit/utils/utils.py @@ -345,17 +345,16 @@ def fluxesFromFootprints(footprints, parentImage, subtractImageMedian=False): badTypeMsg = ("This function works with FootprintSets, single Footprints, and iterables of Footprints. " f"Got {type(footprints)}: {footprints}") if isinstance(footprints, FootprintSet): - fps = footprints.getFootprints() + footprints = footprints.getFootprints() elif isinstance(footprints, Iterable): if not isinstance(footprints[0], Footprint): raise TypeError(badTypeMsg) - fps = footprints elif isinstance(footprints, Footprint): - fps = [footprints] + footprints = [footprints] else: raise TypeError(badTypeMsg) - return [fluxFromFootprint(fp, parentImage, backgroundValue=median) for fp in fps] + return [fluxFromFootprint(fp, parentImage, backgroundValue=median) for fp in footprints] def fluxFromFootprint(footprint, parentImage, backgroundValue=0): From fafac613e2ca3477124d11afa1e6fec9639a7d38 Mon Sep 17 00:00:00 2001 From: Merlin Fisher-Levine Date: Tue, 25 Oct 2022 08:00:48 -0700 Subject: [PATCH 5/7] Change from allclose to all identical --- tests/test_utils.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_utils.py b/tests/test_utils.py index 98e7c737..076e141e 100644 --- a/tests/test_utils.py +++ b/tests/test_utils.py @@ -251,7 +251,7 @@ def test_fluxFromFootprint(self): # functionally unchanged oldExpArray = copy.deepcopy(self.exp.image.array) fluxes = fluxesFromFootprints(obj, self.exp, subtractImageMedian=True) - self.assertTrue(np.allclose(self.exp.image.array, oldExpArray)) + self.assertTrue(np.alltrue(np.equal(self.exp.image.array, oldExpArray))) class TestMemory(lsst.utils.tests.MemoryTestCase): From 4ec63cf890b196d82c37577fafcb8290993100af Mon Sep 17 00:00:00 2001 From: Merlin Fisher-Levine Date: Tue, 25 Oct 2022 09:40:11 -0700 Subject: [PATCH 6/7] Pull tests out to be individual --- tests/test_utils.py | 30 +++++++++++++++++++++--------- 1 file changed, 21 insertions(+), 9 deletions(-) diff --git a/tests/test_utils.py b/tests/test_utils.py index 076e141e..b1b80dbb 100644 --- a/tests/test_utils.py +++ b/tests/test_utils.py @@ -50,7 +50,6 @@ from lsst.summit.utils.butlerUtils import makeDefaultLatissButler from lsst.obs.lsst.translators.latiss import AUXTEL_LOCATION -from lsst.afw.detection import FootprintSet class ExpSkyPositionOffsetTestCase(lsst.utils.tests.TestCase): @@ -238,19 +237,32 @@ def test_fluxFromFootprint(self): # check it can accept a footprintSet, and single and iterables of # footprints - for obj in (footPrintSet, singleFootprint, twoFootprints): - fluxes = fluxesFromFootprints(obj, self.exp) - if isinstance(obj, FootprintSet): - expectedLength = len(footPrintSet.getFootprints()) - self.assertEqual(len(fluxes), expectedLength) # always one flux per footprint - self.assertIsInstance(fluxes, Iterable) - self.assertIsInstance(fluxes[0], np.floating) # always returns floats + # check the footPrintSet + fluxes = fluxesFromFootprints(footPrintSet, self.exp) + expectedLength = len(footPrintSet.getFootprints()) + self.assertEqual(len(fluxes), expectedLength) # always one flux per footprint + self.assertIsInstance(fluxes, Iterable) + self.assertIsInstance(fluxes[0], np.floating) # always returns floats + + # check the singleFootprint + fluxes = fluxesFromFootprints(singleFootprint, self.exp) + expectedLength = 1 + self.assertEqual(len(fluxes), expectedLength) # always one flux per footprint + self.assertIsInstance(fluxes, Iterable) + self.assertIsInstance(fluxes[0], np.floating) # always returns floats + + # check the list of twoFootprints + fluxes = fluxesFromFootprints(twoFootprints, self.exp) + expectedLength = 2 + self.assertEqual(len(fluxes), expectedLength) # always one flux per footprint + self.assertIsInstance(fluxes, Iterable) + self.assertIsInstance(fluxes[0], np.floating) # always returns floats # ensure that subtracting the image median from fluxes leave image # functionally unchanged oldExpArray = copy.deepcopy(self.exp.image.array) - fluxes = fluxesFromFootprints(obj, self.exp, subtractImageMedian=True) + fluxes = fluxesFromFootprints(footPrintSet, self.exp, subtractImageMedian=True) self.assertTrue(np.alltrue(np.equal(self.exp.image.array, oldExpArray))) From 287fe204faeef6be0c9fffbe31b5adfddb93ba51 Mon Sep 17 00:00:00 2001 From: fred3m Date: Mon, 22 May 2023 10:39:26 -0700 Subject: [PATCH 7/7] Update with new afw footprint API --- python/lsst/summit/utils/utils.py | 17 +++---- tests/test_utils.py | 75 ++++++++++++++++--------------- 2 files changed, 48 insertions(+), 44 deletions(-) diff --git a/python/lsst/summit/utils/utils.py b/python/lsst/summit/utils/utils.py index 00440de4..3088ba1a 100644 --- a/python/lsst/summit/utils/utils.py +++ b/python/lsst/summit/utils/utils.py @@ -321,8 +321,8 @@ def fluxesFromFootprints(footprints, parentImage, subtractImageMedian=False): `lsst.afw.detection.Footprint` or `iterable` of `lsst.afw.detection.Footprint` The footprints to measure. - parentImage : `lsst.afw.image.Exposure` - The parent exposure. + parentImage : `lsst.afw.image.Image` + The parent image. subtractImageMedian : `bool`, optional Subtract a whole-image median from each pixel in the footprint when summing as a very crude background subtraction. Does not change the @@ -339,7 +339,7 @@ def fluxesFromFootprints(footprints, parentImage, subtractImageMedian=False): """ median = 0 if subtractImageMedian: - median = np.nanmedian(parentImage.image.array) + median = np.nanmedian(parentImage.array) # poor person's single dispatch badTypeMsg = ("This function works with FootprintSets, single Footprints, and iterables of Footprints. " @@ -354,7 +354,7 @@ def fluxesFromFootprints(footprints, parentImage, subtractImageMedian=False): else: raise TypeError(badTypeMsg) - return [fluxFromFootprint(fp, parentImage, backgroundValue=median) for fp in footprints] + return np.array([fluxFromFootprint(fp, parentImage, backgroundValue=median) for fp in footprints]) def fluxFromFootprint(footprint, parentImage, backgroundValue=0): @@ -366,8 +366,8 @@ def fluxFromFootprint(footprint, parentImage, backgroundValue=0): ---------- footprint : `lsst.afw.detection.Footprint` The footprint to measure. - parentImage : `lsst.afw.image.Exposure` - The parent exposure. + parentImage : `lsst.afw.image.Image` + Image containing the footprint. backgroundValue : `bool`, optional The value to subtract from each pixel in the footprint when summing as a very crude background subtraction. Does not change the original @@ -379,8 +379,9 @@ def fluxFromFootprint(footprint, parentImage, backgroundValue=0): The flux in the footprint """ if backgroundValue: # only do the subtraction if non-zero for speed - return footprint.getSpans().flatten(parentImage.image.array - backgroundValue).sum() - return footprint.getSpans().flatten(parentImage.image.array).sum() + xy0 = parentImage.getBBox().getMin() + return footprint.computeFluxFromArray(parentImage.array - backgroundValue, xy0) + return footprint.computeFluxFromImage(parentImage) def humanNameForCelestialObject(objName): diff --git a/tests/test_utils.py b/tests/test_utils.py index b1b80dbb..4715bb6e 100644 --- a/tests/test_utils.py +++ b/tests/test_utils.py @@ -28,7 +28,9 @@ import astropy.time import astropy.units as u +import lsst.afw.detection as afwDetect import lsst.afw.image as afwImage +import lsst.afw.geom as afwGeom import lsst.geom as geom import lsst.utils.tests import numpy as np @@ -44,10 +46,7 @@ quickSmooth, getQuantiles, fluxesFromFootprints, - detectObjectsInExp, ) -from lsst.summit.utils.bestEffort import BestEffortIsr -from lsst.summit.utils.butlerUtils import makeDefaultLatissButler from lsst.obs.lsst.translators.latiss import AUXTEL_LOCATION @@ -215,55 +214,59 @@ def test_quantiles(self): class ImageBasedTestCase(lsst.utils.tests.TestCase): - """A test case for testing sky position offsets for exposures.""" + def test_fluxFromFootprint(self): + image = afwImage.Image( + np.arange(8100, dtype=np.int32).reshape(90, 90), + xy0=lsst.geom.Point2I(10, 12), + dtype="I" + ) - @classmethod - def setUpClass(cls): - try: - cls.butler = makeDefaultLatissButler() - except FileNotFoundError: - raise unittest.SkipTest("Skipping tests that require the LATISS butler repo.") + radius = 3 + spans = afwGeom.SpanSet.fromShape(radius, afwGeom.Stencil.CIRCLE, offset=(27, 30)) + footprint1 = afwDetect.Footprint(spans) - cls.dataId = {'day_obs': 20200315, 'seq_num': 120, 'detector': 0} - cls.bestEffort = BestEffortIsr() # should always succeed if makeDefaultLatissButler works - cls.exp = cls.bestEffort.getExposure(cls.dataId) + # The extracted footprint should be the same as the product of the + # spans and the overlapped bow with the image + truth1 = spans.asArray() * image.array[15:22, 14:21] - def test_fluxFromFootprint(self): - footPrintSet = detectObjectsInExp(self.exp) - allFootprints = footPrintSet.getFootprints() - self.assertGreaterEqual(len(allFootprints), 3) - singleFootprint = allFootprints[0] - twoFootprints = allFootprints[0:2] + radius = 3 + spans = afwGeom.SpanSet.fromShape(radius, afwGeom.Stencil.CIRCLE, offset=(44, 49)) + footprint2 = afwDetect.Footprint(spans) + truth2 = spans.asArray() * image.array[34:41, 31:38] + + allFootprints = [footprint1, footprint2] + footprintSet = afwDetect.FootprintSet(image.getBBox()) + footprintSet.setFootprints(allFootprints) # check it can accept a footprintSet, and single and iterables of # footprints + with self.assertRaises(TypeError): + fluxesFromFootprints(10, image) - # check the footPrintSet - fluxes = fluxesFromFootprints(footPrintSet, self.exp) - expectedLength = len(footPrintSet.getFootprints()) - self.assertEqual(len(fluxes), expectedLength) # always one flux per footprint - self.assertIsInstance(fluxes, Iterable) - self.assertIsInstance(fluxes[0], np.floating) # always returns floats + with self.assertRaises(TypeError): + fluxesFromFootprints([8, 6, 7, 5, 3, 0, 9], image) - # check the singleFootprint - fluxes = fluxesFromFootprints(singleFootprint, self.exp) - expectedLength = 1 + # check the footPrintSet + fluxes = fluxesFromFootprints(footprintSet, image) + expectedLength = len(footprintSet.getFootprints()) self.assertEqual(len(fluxes), expectedLength) # always one flux per footprint self.assertIsInstance(fluxes, Iterable) - self.assertIsInstance(fluxes[0], np.floating) # always returns floats + self.assertAlmostEqual(fluxes[0], np.sum(truth1)) + self.assertAlmostEqual(fluxes[1], np.sum(truth2)) - # check the list of twoFootprints - fluxes = fluxesFromFootprints(twoFootprints, self.exp) + # check the list of footprints + fluxes = fluxesFromFootprints(allFootprints, image) expectedLength = 2 self.assertEqual(len(fluxes), expectedLength) # always one flux per footprint self.assertIsInstance(fluxes, Iterable) - self.assertIsInstance(fluxes[0], np.floating) # always returns floats + self.assertAlmostEqual(fluxes[0], np.sum(truth1)) + self.assertAlmostEqual(fluxes[1], np.sum(truth2)) # ensure that subtracting the image median from fluxes leave image - # functionally unchanged - oldExpArray = copy.deepcopy(self.exp.image.array) - fluxes = fluxesFromFootprints(footPrintSet, self.exp, subtractImageMedian=True) - self.assertTrue(np.alltrue(np.equal(self.exp.image.array, oldExpArray))) + # pixels untouched + oldImageArray = copy.deepcopy(image.array) + fluxes = fluxesFromFootprints(footprintSet, image, subtractImageMedian=True) + np.testing.assert_array_equal(image.array, oldImageArray) class TestMemory(lsst.utils.tests.MemoryTestCase):