Skip to content

Commit

Permalink
Adding matching kernel spatially sampled metrics
Browse files Browse the repository at this point in the history
  • Loading branch information
@BrunoSanchez
BrunoSanchez committed Jun 19, 2024
1 parent 7202c41 commit 17542f2
Show file tree
Hide file tree
Showing 2 changed files with 72 additions and 16 deletions.
63 changes: 59 additions & 4 deletions python/lsst/ip/diffim/computeSpatiallySampledMetrics.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,7 @@
import lsst.pipe.base as pipeBase
import lsst.pex.config as pexConfig

from lsst.ip.diffim.utils import getPsfFwhm, angleMean, evaluateMaskFraction
from lsst.ip.diffim.utils import getPsfFwhm, angleMean, evaluateMaskFraction, getKernelCenterDisplacement
from lsst.meas.algorithms import SkyObjectsTask
from lsst.pex.exceptions import InvalidParameterError, RangeError
from lsst.utils.timer import timeMethod
Expand Down Expand Up @@ -72,6 +72,12 @@ class SpatiallySampledMetricsConnections(pipeBase.PipelineTaskConnections,
storageClass="SourceCatalog",
name="{fakesType}{coaddName}Diff_candidateDiaSrc",
)
psfMatchingKernel = pipeBase.connectionTypes.Input(
doc="Kernel used to PSF match the science and template images.",
dimensions=("instrument", "visit", "detector"),
storageClass="MatchingKernel",
name="{fakesType}{coaddName}Diff_psfMatchKernel",
)
spatiallySampledMetrics = pipeBase.connectionTypes.Output(
doc="Summary metrics computed at randomized locations.",
dimensions=("instrument", "visit", "detector"),
Expand Down Expand Up @@ -164,9 +170,32 @@ def __init__(self, **kwargs):
"%s_mask_fraction"%maskPlane.lower(), "F",
"Fraction of pixels with %s mask"%maskPlane
)
self.schema.addField(
"psfMatchingKernel_sum", "F",
"PSF matching kernel sum at location.")
self.schema.addField(
"psfMatchingKernel_dx", "F",
"PSF matching kernel centroid offset in x at location.",
units="pixel")
self.schema.addField(
"psfMatchingKernel_dy", "F",
"PSF matching kernel centroid offset in y at location.",
units="pixel")
self.schema.addField(
"psfMatchingKernel_length", "F",
"PSF matching kernel centroid offset module.",
units="arcsecond")
self.schema.addField(
"psfMatchingKernel_position_angle", "F",
"PSF matching kernel centroid offset position angle.",
units="radian")
self.schema.addField(
"psfMatchingKernel_direction", "F",
"PSF matching kernel centroid offset direction in detector plane.",
units="radian")

@timeMethod
def run(self, science, matchedTemplate, template, difference, diaSources):
def run(self, science, matchedTemplate, template, difference, diaSources, psfMatchingKernel):
"""Calculate difference image metrics on specific locations across the images
Parameters
Expand All @@ -183,6 +212,8 @@ def run(self, science, matchedTemplate, template, difference, diaSources):
Result of subtracting template from the science image.
diaSources : `lsst.afw.table.SourceCatalog`
The catalog of detected sources.
psfMatchingKernel : `~lsst.afw.math.LinearCombinationKernel`
The PSF matching kernel of the subtraction to evaluate.
Returns
-------
Expand All @@ -207,12 +238,13 @@ def run(self, science, matchedTemplate, template, difference, diaSources):

for src in spatiallySampledMetrics:
self._evaluateLocalMetric(src, science, matchedTemplate, template, difference, diaSources,
metricsMaskPlanes=metricsMaskPlanes)
metricsMaskPlanes=metricsMaskPlanes,
psfMatchingKernel=psfMatchingKernel)

return pipeBase.Struct(spatiallySampledMetrics=spatiallySampledMetrics.asAstropy())

def _evaluateLocalMetric(self, src, science, matchedTemplate, template, difference, diaSources,
metricsMaskPlanes):
metricsMaskPlanes, psfMatchingKernel):
"""Calculate image quality metrics at spatially sampled locations.
Parameters
Expand All @@ -229,6 +261,8 @@ def _evaluateLocalMetric(self, src, science, matchedTemplate, template, differen
Result of subtracting template from the science image.
metricsMaskPlanes : `list` of `str`
Mask planes to calculate metrics from.
psfMatchingKernel : `~lsst.afw.math.LinearCombinationKernel`
The PSF matching kernel of the subtraction to evaluate.
"""
bbox = src.getFootprint().getBBox()
pix = bbox.getCenter()
Expand Down Expand Up @@ -271,3 +305,24 @@ def _evaluateLocalMetric(self, src, science, matchedTemplate, template, differen
src.set("%s_mask_fraction"%maskPlane.lower(),
evaluateMaskFraction(difference.mask[bbox], maskPlane)
)

krnlSum, dx, dy, direction, length = getKernelCenterDisplacement(
psfMatchingKernel, src.get('x'), src.get('y'))

point1 = lsst.geom.SpherePoint(
src.get('coord_ra'), src.get('coord_dec'),
lsst.geom.radians)
point2 = science.wcs.pixelToSky(src.get('x') + dx, src.get('y') + dy)
bearing = point1.bearingTo(point2)
pa_ref_angle = lsst.geom.Angle(np.pi/2, lsst.geom.radians)
pa = pa_ref_angle - bearing
# Wrap around to get Delta_RA from -pi to +pi
pa = pa.wrapCtr()
position_angle = pa.asRadians()

src.set('psfMatchingKernel_sum', krnlSum)
src.set('psfMatchingKernel_dx', dx)
src.set('psfMatchingKernel_dy', dy)
src.set('psfMatchingKernel_length', length*pixScale.asArcseconds())
src.set('psfMatchingKernel_position_angle', position_angle) # in E of N position angle
src.set('psfMatchingKernel_direction', direction) # direction offset in detector
25 changes: 13 additions & 12 deletions python/lsst/ip/diffim/utils.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,7 @@
"""Support utilities for Measuring sources"""

# Export DipoleTestImage to expose fake image generating funcs
__all__ = ["DipoleTestImage", "evaluateMeanPsfFwhm", "getPsfFwhm"]
__all__ = ["DipoleTestImage", "evaluateMeanPsfFwhm", "getPsfFwhm", "getKernelCenterDisplacement"]

import itertools
import numpy as np
Expand Down Expand Up @@ -853,12 +853,9 @@ def detectDipoleSources(self, doMerge=True, diffim=None, detectSigma=5.5, grow=3
return detectTask, schema


def getKernelCenterDisplacement(kernel, x, y, im=None):
def getKernelCenterDisplacement(kernel, x, y, image=None):
"""Calculate the PSF matching kernel peak offset from the nominal
position.
This is following Robert's code at
https://lsstc.slack.com/archives/C2JPMCF5X/p1715784207956889
https://lsstc.slack.com/archives/C2JPMCF5X/p1715865918632319
Parameters
----------
Expand All @@ -868,29 +865,33 @@ def getKernelCenterDisplacement(kernel, x, y, im=None):
The x position on the detector to evaluate the kernel
y : `float`
The y position on the detector to evaluate the kernel
im : `~lsst.afw.image._image.ImageD`
image : `~lsst.afw.image.ImageD`
The image to use as base for computing kernel pixel values
Returns
-------
krnl_sum : `float`
kernel_sum : `float`
The sum of the kernel on the desired location
dx : `float`
The displacement of the kernel averaged peak, with respect to the
center of the extraction of the kernel
dy : `float`
The displacement of the kernel averaged peak, with respect to the
center of the extraction of the kernel
pos_angle: `float`
The position angle in detector coordinates of the displacement
length : `float`
The displacement module of the kernel centroid in pixel units
"""

if im is None:
im = afwImage.ImageD(kernel.getDimensions())
if image is None:
image = afwImage.ImageD(kernel.getDimensions())

# obtain the kernel image
hsize = kernel.getWidth()//2
krnl_sum = kernel.computeImage(im, doNormalize=False, x=x, y=y)
kernel_sum = kernel.computeImage(image, doNormalize=False, x=x, y=y)

data = im.array
data = image.array
h, w = data.shape
xx = np.arange(w)
yy = np.arange(h)
Expand All @@ -908,7 +909,7 @@ def getKernelCenterDisplacement(kernel, x, y, im=None):
pos_angle = np.arctan2(dy, dx)
length = np.sqrt(dx**2 + dy**2)

return krnl_sum, dx, dy, pos_angle, length
return kernel_sum, dx, dy, pos_angle, length


def getPsfFwhm(psf, average=True, position=None):
Expand Down

0 comments on commit 17542f2

Please sign in to comment.