From 1c9dda1e8135f2fbe290ec1755abb6e587263d2d Mon Sep 17 00:00:00 2001 From: Johnny Esteves Date: Wed, 25 Feb 2026 10:30:29 -0800 Subject: [PATCH 1/6] Fix FootprintSet STDEV detection failure in runSourceDetection Replace detectObjectsInExp with direct detection using Threshold.VALUE and np.std() to avoid STDEVCLIP returning 0 for some guider images. Add starMaskThreshold parameter to exclude bright pixels when computing background standard deviation. --- python/lsst/summit/utils/guiders/detection.py | 34 +++++++++++++++---- 1 file changed, 28 insertions(+), 6 deletions(-) diff --git a/python/lsst/summit/utils/guiders/detection.py b/python/lsst/summit/utils/guiders/detection.py index 4f17b3c2..e109398c 100644 --- a/python/lsst/summit/utils/guiders/detection.py +++ b/python/lsst/summit/utils/guiders/detection.py @@ -37,8 +37,8 @@ from astropy.nddata import Cutout2D from astropy.stats import sigma_clipped_stats +import lsst.afw.detection as afwDetect from lsst.afw.image import ExposureF, ImageF, MaskedImageF -from lsst.summit.utils.utils import detectObjectsInExp from .reading import GuiderData @@ -277,6 +277,7 @@ def runSourceDetection( cutOutSize: int = 25, apertureRadius: int = 5, gain: float = 1.0, + starMaskThreshold: float = 20000.0, ) -> pd.DataFrame: """ Detect sources in an image and measure their properties. @@ -293,6 +294,9 @@ def runSourceDetection( Aperture radius in pixels for photometry. gain : `float` Detector gain (e-/ADU). + starMaskThreshold : `float` + Pixel count threshold for masking bright stars when computing the + background standard deviation. Pixels above this value are excluded. Returns ------- @@ -303,12 +307,30 @@ def runSourceDetection( exposure = ExposureF(MaskedImageF(ImageF(image))) # Step 2: Detect sources - # we assume that we have bright stars - # filter out stamps with no stars + # We use Threshold.VALUE with np.std() instead of Threshold.STDEV because + # STDEV uses STDEVCLIP internally which can return 0 for some guider + # images, causing FootprintSet to fail with "St. dev. must be > 0". + footprints = None if not isBlankImage(image): - footprints = detectObjectsInExp(exposure, nSigma=threshold) - else: - footprints = None + median = np.nanmedian(exposure.image.array) + exposure.image -= median + # Mask bright stars when computing std to get background noise estimate + bgPixels = exposure.image.array[exposure.image.array < starMaskThreshold] + if len(bgPixels) == 0: + raise ValueError( + f"No pixels below starMaskThreshold={starMaskThreshold}. " + "All pixels are saturated or threshold is too low." + ) + imageStd = np.std(bgPixels) + if imageStd <= 0: + raise ValueError( + f"Background standard deviation is {imageStd}. " + "Image may be flat or starMaskThreshold may need adjustment." + ) + absThreshold = threshold * imageStd + thresh = afwDetect.Threshold(absThreshold, afwDetect.Threshold.VALUE) + footprints = afwDetect.FootprintSet(exposure.getMaskedImage(), thresh, "DETECTED", 10) + exposure.image += median if not footprints: return pd.DataFrame(columns=DEFAULT_COLUMNS) From c0414f51966e56716929e12aa6fa603fec16e1ad Mon Sep 17 00:00:00 2001 From: Johnny Esteves Date: Mon, 2 Mar 2026 16:17:51 -0800 Subject: [PATCH 2/6] Suppress pre-existing mypy union-attr errors in utils.py --- python/lsst/summit/utils/utils.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/python/lsst/summit/utils/utils.py b/python/lsst/summit/utils/utils.py index 471ef6ea..bf3a2fe4 100644 --- a/python/lsst/summit/utils/utils.py +++ b/python/lsst/summit/utils/utils.py @@ -488,8 +488,8 @@ def _getAltAzZenithsFromSeqNum( for seqNum in seqNumList: md = butler.get("raw.metadata", day_obs=dayObs, seq_num=seqNum, detector=0) obsInfo = ObservationInfo(md) - alt = obsInfo.altaz_begin.alt.value - az = obsInfo.altaz_begin.az.value + alt = obsInfo.altaz_begin.alt.value # type: ignore[union-attr] + az = obsInfo.altaz_begin.az.value # type: ignore[union-attr] elevations.append(alt) zeniths.append(90 - alt) azimuths.append(az) From cce93dc37926b73186fd831ac6b496a03be183f3 Mon Sep 17 00:00:00 2001 From: Johnny Esteves Date: Mon, 2 Mar 2026 16:25:06 -0800 Subject: [PATCH 3/6] Add nPixMin parameter to runSourceDetection --- python/lsst/summit/utils/guiders/detection.py | 34 +++++++------------ 1 file changed, 13 insertions(+), 21 deletions(-) diff --git a/python/lsst/summit/utils/guiders/detection.py b/python/lsst/summit/utils/guiders/detection.py index e109398c..b399e9ce 100644 --- a/python/lsst/summit/utils/guiders/detection.py +++ b/python/lsst/summit/utils/guiders/detection.py @@ -277,7 +277,7 @@ def runSourceDetection( cutOutSize: int = 25, apertureRadius: int = 5, gain: float = 1.0, - starMaskThreshold: float = 20000.0, + nPixMin: int = 10, ) -> pd.DataFrame: """ Detect sources in an image and measure their properties. @@ -294,9 +294,8 @@ def runSourceDetection( Aperture radius in pixels for photometry. gain : `float` Detector gain (e-/ADU). - starMaskThreshold : `float` - Pixel count threshold for masking bright stars when computing the - background standard deviation. Pixels above this value are excluded. + nPixMin : `int` + Minimum number of pixels in a footprint for detection. Returns ------- @@ -307,29 +306,22 @@ def runSourceDetection( exposure = ExposureF(MaskedImageF(ImageF(image))) # Step 2: Detect sources - # We use Threshold.VALUE with np.std() instead of Threshold.STDEV because - # STDEV uses STDEVCLIP internally which can return 0 for some guider + # We use Threshold.VALUE with a sigma68 background estimate instead of + # Threshold.STDEV because STDEVCLIP can return 0 for some guider # images, causing FootprintSet to fail with "St. dev. must be > 0". + # sigma68 = (p84 - p16) / 2 is robust to bright stars (they only + # affect the upper tail) and doesn't fail on flat-background images + # the way MAD does (MAD=0 when >50% of pixels are identical). footprints = None if not isBlankImage(image): - median = np.nanmedian(exposure.image.array) - exposure.image -= median - # Mask bright stars when computing std to get background noise estimate - bgPixels = exposure.image.array[exposure.image.array < starMaskThreshold] - if len(bgPixels) == 0: - raise ValueError( - f"No pixels below starMaskThreshold={starMaskThreshold}. " - "All pixels are saturated or threshold is too low." - ) - imageStd = np.std(bgPixels) + p16, median, p84 = np.nanpercentile(image, [16, 50, 84]) + imageStd = (p84 - p16) / 2.0 if imageStd <= 0: - raise ValueError( - f"Background standard deviation is {imageStd}. " - "Image may be flat or starMaskThreshold may need adjustment." - ) + return pd.DataFrame(columns=DEFAULT_COLUMNS) + exposure.image -= median absThreshold = threshold * imageStd thresh = afwDetect.Threshold(absThreshold, afwDetect.Threshold.VALUE) - footprints = afwDetect.FootprintSet(exposure.getMaskedImage(), thresh, "DETECTED", 10) + footprints = afwDetect.FootprintSet(exposure.getMaskedImage(), thresh, "DETECTED", nPixMin) exposure.image += median if not footprints: From a887584fef09be9dea1c81de0ec1ab05ae833990 Mon Sep 17 00:00:00 2001 From: Johnny Esteves Date: Mon, 2 Mar 2026 18:13:53 -0800 Subject: [PATCH 4/6] Add dither to coadd median to prevent STDEVCLIP=0 Add a uniform [0, 1) random dither per stamp before computing the median coadd. This breaks integer quantization that causes the pixel distribution to be so narrow that STDEVCLIP returns 0 on some guider images. Fix suggested by Robert Lupton. (DM-54263) --- python/lsst/summit/utils/guiders/reading.py | 11 ++++++++--- 1 file changed, 8 insertions(+), 3 deletions(-) diff --git a/python/lsst/summit/utils/guiders/reading.py b/python/lsst/summit/utils/guiders/reading.py index d671fb5f..2ca65206 100644 --- a/python/lsst/summit/utils/guiders/reading.py +++ b/python/lsst/summit/utils/guiders/reading.py @@ -406,10 +406,15 @@ def getStampArrayCoadd(self, detName: str) -> np.ndarray: if len(stamps) == 0: raise ValueError(f"No stamps found for detector {detName!r}") - # Collect arrays, with optional bias subtraction + # Collect arrays, with optional bias subtraction. arrList = [self[detName, idx] for idx in range(len(stamps))] - stack = np.nanmedian(arrList, axis=0) - return stack + # Add a uniform [0, 1) dither per stamp before the median to break + # integer quantization. Without this, the coadd pixel distribution + # can be so narrow that STDEVCLIP returns 0. (See DM-54263.) + stack = np.array(arrList, dtype=np.float32) + rng = np.random.default_rng(seed=0) + stack += rng.uniform(0, 1, size=stack.shape).astype(np.float32) + return np.nanmedian(stack, axis=0) def getGuiderAmpName(self, detName: str) -> str: """ From 0a211d2dfcfa1f8f3356ff2635dd3e4c6844185d Mon Sep 17 00:00:00 2001 From: Johnny Esteves Date: Mon, 2 Mar 2026 19:08:57 -0800 Subject: [PATCH 5/6] Use STDEVCLIP with sigma68 fallback for source detection Replace the manual sigma68-only background estimation with STDEVCLIP from lsst.afw.math, falling back to sigma68 if STDEVCLIP returns 0. The dithered coadds (previous commit) prevent STDEVCLIP from failing, but the fallback provides a safety net. Also adds a configurable nPixMin parameter for footprint detection. (DM-54263) --- python/lsst/summit/utils/guiders/detection.py | 24 ++++++++++--------- 1 file changed, 13 insertions(+), 11 deletions(-) diff --git a/python/lsst/summit/utils/guiders/detection.py b/python/lsst/summit/utils/guiders/detection.py index b399e9ce..a2c95bdc 100644 --- a/python/lsst/summit/utils/guiders/detection.py +++ b/python/lsst/summit/utils/guiders/detection.py @@ -39,6 +39,7 @@ import lsst.afw.detection as afwDetect from lsst.afw.image import ExposureF, ImageF, MaskedImageF +from lsst.afw.math import STDEVCLIP, makeStatistics from .reading import GuiderData @@ -305,20 +306,22 @@ def runSourceDetection( # Step 1: Convert numpy image to MaskedImage and Exposure exposure = ExposureF(MaskedImageF(ImageF(image))) - # Step 2: Detect sources - # We use Threshold.VALUE with a sigma68 background estimate instead of - # Threshold.STDEV because STDEVCLIP can return 0 for some guider - # images, causing FootprintSet to fail with "St. dev. must be > 0". - # sigma68 = (p84 - p16) / 2 is robust to bright stars (they only - # affect the upper tail) and doesn't fail on flat-background images - # the way MAD does (MAD=0 when >50% of pixels are identical). + # Step 2: Detect sources using STDEVCLIP for the background noise. + # The input coadd images are dithered (see GuiderData.getStampArrayCoadd) + # to prevent integer quantization from collapsing the pixel distribution + # and causing STDEVCLIP to return 0. (See DM-54263.) footprints = None if not isBlankImage(image): - p16, median, p84 = np.nanpercentile(image, [16, 50, 84]) - imageStd = (p84 - p16) / 2.0 + median = np.nanmedian(image) + exposure.image -= median + imageStd = float(makeStatistics(exposure.getMaskedImage(), STDEVCLIP).getValue(STDEVCLIP)) + if imageStd <= 0: + # Fallback: sigma68 is robust to quantization and bright stars. + p16, p84 = np.nanpercentile(image, [16, 84]) + imageStd = (p84 - p16) / 2.0 if imageStd <= 0: + exposure.image += median return pd.DataFrame(columns=DEFAULT_COLUMNS) - exposure.image -= median absThreshold = threshold * imageStd thresh = afwDetect.Threshold(absThreshold, afwDetect.Threshold.VALUE) footprints = afwDetect.FootprintSet(exposure.getMaskedImage(), thresh, "DETECTED", nPixMin) @@ -627,7 +630,6 @@ def buildReferenceCatalog( apertureRadius = int(config.aperSizeArcsec / pixelScale) array = guiderData.getStampArrayCoadd(guiderName) - # array = np.where(array < 0, 0, array) # Ensure no negative values array = array - np.nanmin(array) # Ensure no negative values sources = runSourceDetection( array, From dc35eb5128e13116ecd745724659e4c12df543c3 Mon Sep 17 00:00:00 2001 From: Johnny Esteves Date: Tue, 3 Mar 2026 06:11:28 -0800 Subject: [PATCH 6/6] Replace type: ignore with assertion for altaz_begin Address review feedback from Merlin: use explicit assertion instead of mypy ignore comment. (DM-54263) --- python/lsst/summit/utils/utils.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/python/lsst/summit/utils/utils.py b/python/lsst/summit/utils/utils.py index bf3a2fe4..ec8750e1 100644 --- a/python/lsst/summit/utils/utils.py +++ b/python/lsst/summit/utils/utils.py @@ -488,8 +488,9 @@ def _getAltAzZenithsFromSeqNum( for seqNum in seqNumList: md = butler.get("raw.metadata", day_obs=dayObs, seq_num=seqNum, detector=0) obsInfo = ObservationInfo(md) - alt = obsInfo.altaz_begin.alt.value # type: ignore[union-attr] - az = obsInfo.altaz_begin.az.value # type: ignore[union-attr] + assert obsInfo.altaz_begin is not None, f"altaz_begin is None for seqNum={seqNum}" + alt = obsInfo.altaz_begin.alt.value + az = obsInfo.altaz_begin.az.value elevations.append(alt) zeniths.append(90 - alt) azimuths.append(az)