From b5a2b3996d3fd1a0eb4714dac6c5df9906d82ebf Mon Sep 17 00:00:00 2001 From: Clare Saunders Date: Thu, 29 Jan 2026 07:51:01 -0800 Subject: [PATCH] Calculate image/patch overlap in patch coordinates --- python/lsst/ip/diffim/getTemplate.py | 13 +++++++++---- 1 file changed, 9 insertions(+), 4 deletions(-) diff --git a/python/lsst/ip/diffim/getTemplate.py b/python/lsst/ip/diffim/getTemplate.py index 70382b038..6d03a25bf 100644 --- a/python/lsst/ip/diffim/getTemplate.py +++ b/python/lsst/ip/diffim/getTemplate.py @@ -243,6 +243,7 @@ def getExposures(self, coaddExposureHandles, bbox, skymap, wcs): # Exposure's validPolygon would be more accurate detectorPolygon = geom.Box2D(bbox) + detectorCorners = wcs.pixelToSky(detectorPolygon.getCorners()) overlappingArea = 0 coaddExposures = collections.defaultdict(list) dataIds = collections.defaultdict(list) @@ -251,11 +252,15 @@ def getExposures(self, coaddExposureHandles, bbox, skymap, wcs): dataId = coaddRef.dataId patchWcs = skymap[dataId["tract"]].getWcs() patchBBox = skymap[dataId["tract"]][dataId["patch"]].getOuterBBox() - patchCorners = patchWcs.pixelToSky(geom.Box2D(patchBBox).getCorners()) - patchPolygon = afwGeom.Polygon(wcs.skyToPixel(patchCorners)) - if patchPolygon.intersection(detectorPolygon): + patchPolygon = afwGeom.Polygon(geom.Box2D(patchBBox)) + # Calculate detector/patch overlap in patch coordinates rather than + # detector coordinates because the skymap's inverse mapping + # (patchWcs.skyToPixel()) is more stable than the detector's for + # arbitrary sky coordinates. + detectorInPatchCoordinates = afwGeom.Polygon(patchWcs.skyToPixel(detectorCorners)) + if patchPolygon.intersection(detectorInPatchCoordinates): overlappingArea += patchPolygon.intersectionSingle( - detectorPolygon + detectorInPatchCoordinates ).calculateArea() self.log.info( "Using template input tract=%s, patch=%s",