Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

DM-29673 Add centre of mass centroiding to support donut images #497

Merged
merged 6 commits into from
Apr 14, 2021
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
120 changes: 84 additions & 36 deletions python/lsst/pipe/tasks/quickFrameMeasurement.py
Original file line number Diff line number Diff line change
Expand Up @@ -20,8 +20,11 @@
# along with this program. If not, see <https://www.gnu.org/licenses/>.

import numpy as np
import scipy.ndimage as ndImage

import lsst.afw.detection as afwDetect
import lsst.afw.table as afwTable
import lsst.geom as geom
import lsst.meas.base as measBase
import lsst.daf.base as dafBase
import lsst.pipe.base as pipeBase
Expand Down Expand Up @@ -79,6 +82,11 @@ class QuickFrameMeasurementTaskConfig(pexConfig.Config):
doc="Minimum number of pixels in a detected source.",
default=10,
)
donutDiameter = pexConfig.Field(
dtype=int,
doc="The expected diameter of donuts in a donut image, in pixels.",
default=400,
)

def setDefaults(self):
super().setDefaults()
Expand Down Expand Up @@ -107,15 +115,14 @@ class QuickFrameMeasurementTask(pipeBase.Task):

Parameters
----------
config : lsst.pipe.tasks.quickFrameMeasurement.QuickFrameMeasurementTaskConfig
config : `lsst.pipe.tasks.quickFrameMeasurement.QuickFrameMeasurementTaskConfig`
Configuration class for the QuickFrameMeasurementTask.

display : lsst.afw.display.Display, optional
display : `lsst.afw.display.Display`, optional
The display to use for showing the images, detections and centroids.

Returns
-------
result : lsst.pipe.base.Struct()
result : `lsst.pipe.base.Struct`
Return strucure containing whether the task was successful, the main
source's centroid, its the aperture fluxes, the ixx and iyy of the
source, and the median ixx, iyy of the detections in the exposure.
Expand Down Expand Up @@ -163,21 +170,18 @@ def detectObjectsInExp(exp, nSigma, nPixMin, grow=0):

Parameters
----------
exp : lsst.afw.image.Exposure
exp : `lsst.afw.image.Exposure`
Image in which to detect objects.

nSigma : float
nSigma : `float`
nSigma above image's stddev at which to set the detection threshold.

nPixMin : int
nPixMin : `int`
Minimum number of pixels for detection.

grow : int
grow : `int`
Grow the detected footprint by this many pixels.

Returns
-------
footPrintSet : lsst.afw.detection.FootprintSet
footPrintSet : `lsst.afw.detection.FootprintSet`
FootprintSet containing the detections.
"""
threshold = afwDetect.Threshold(nSigma, afwDetect.Threshold.STDEV)
Expand All @@ -193,15 +197,14 @@ def checkResult(exp, centroid, srcNum, percentile):

Parameters
----------
exp : lsst.afw.image.Exposure
exp : `lsst.afw.image.Exposure`
The exposure on which to operate
centroid : `tuple` of `float`
Location of the centroid in pixel coordinates

scrNum : int
scrNum : `int`
Number of the source in the source catalog. Only used if the check
is failed, for debug purposes.

percentile : float
percentile : `float`
Image's percentile above which the pixel containing the centroid
must be in order to pass the check.

Expand All @@ -227,6 +230,40 @@ def _calcMedianXxYy(objData):
medianYy = np.nanmedian([element['xx'] for element in objData.values()])
return medianXx, medianYy

@staticmethod
def _getCenterOfMass(exp, nominalCentroid, boxSize):
"""Get the centre of mass around a point in the image.

Parameters
----------
exp : `lsst.afw.image.Exposure`
The exposure in question.
nominalCentroid : `tuple` of `float`
Nominal location of the centroid in pixel coordinates.
boxSize : `int`
The size of the box around the nominalCentroid in which to measure
the centre of mass.

Returns
-------
com : `tuple` of `float`
The locaiton of the centre of mass of the brightest source in pixel
coordinates.
"""
centroidPoint = geom.Point2I(nominalCentroid)
extent = geom.Extent2I(1, 1)
bbox = geom.Box2I(centroidPoint, extent)
bbox = bbox.dilatedBy(int(boxSize//2))
mfisherlevine marked this conversation as resolved.
Show resolved Hide resolved
bbox = bbox.clippedTo(exp.getBBox())
data = exp[bbox].image.array
xy0 = exp[bbox].getXY0()

peak = ndImage.center_of_mass(data)
peak = (peak[1], peak[0]) # numpy coords returned
com = geom.Point2D(xy0)
com.shift(geom.Extent2D(*peak))
return (com[0], com[1])

def _calcBrightestObjSrcNum(self, objData):
"""Find the brightest source which passes the cuts among the sources.

Expand All @@ -237,7 +274,7 @@ def _calcBrightestObjSrcNum(self, objData):

Returns
-------
srcNum : int
srcNum : `int`
The source number of the brightest source which passes the cuts.
"""
max70, max70srcNum = -1, -1
Expand Down Expand Up @@ -293,15 +330,14 @@ def _measureFp(self, fp, exp):

Parameters
----------
fp : lsst.afw.detection.Footprint
fp : `lsst.afw.detection.Footprint`
The footprint to measure.

exp : lsst.afw.image.Exposure
exp : `lsst.afw.image.Exposure`
The footprint's parent exposure.

Returns
-------
src : lsst.afw.table.SourceRecord
src : `lsst.afw.table.SourceRecord`
The source record containing the measurements.
"""
src = self.table.makeRecord()
Expand All @@ -316,12 +352,12 @@ def _getDataFromSrcRecord(self, src):

Parameters
----------
src : lsst.afw.table.SourceRecord
src : `lsst.afw.table.SourceRecord`
The source record from which to extract the measurements.

Returns
-------
srcData : lsst.pipe.base.Struct
srcData : `lsst.pipe.base.Struct`
The struct containing the extracted measurements.
"""
pScale = self.plateScale
Expand All @@ -345,14 +381,14 @@ def _getDataFromFootprintOnly(fp, exp):

Parameters
----------
fp : lsst.afw.detection.Footprint
fp : `lsst.afw.detection.Footprint`
The footprint to measure.
exp : lsst.afw.image.Exposure
exp : `lsst.afw.image.Exposure`
The footprint's parent exposure.

Returns
-------
srcData : lsst.pipe.base.Struct
srcData : `lsst.pipe.base.Struct`
The struct containing the extracted measurements.
"""
xx = fp.getShape().getIxx()
Expand All @@ -373,7 +409,7 @@ def _measurementResultToDict(measurementResult):

Parameters
----------
measurementResult : lsst.afw.table.SourceRecord
measurementResult : `lsst.afw.table.SourceRecord`
The source record to convert to a dict.

Returns
Expand All @@ -396,19 +432,20 @@ def _makeEmptyReturnStruct():

Returns
-------
objData : lsst.pipe.base.Struct
objData : `lsst.pipe.base.Struct`
The default template return structure.
"""
result = pipeBase.Struct()
result.success = False
result.brightestObjCentroid = (np.nan, np.nan)
result.brightestObjCentroidCofM = None
result.brightestObj_xXyY = (np.nan, np.nan)
result.brightestObjApFlux70 = np.nan
result.brightestObjApFlux25 = np.nan
result.medianXxYy = (np.nan, np.nan)
return result

def run(self, exp, doDisplay=False):
def run(self, exp, *, donutDiameter=None, doDisplay=False):
mfisherlevine marked this conversation as resolved.
Show resolved Hide resolved
"""Calculate position, flux and shape of the brightest star in an image.

Given an an assembled (and at least minimally ISRed exposure),
Expand All @@ -417,16 +454,19 @@ def run(self, exp, doDisplay=False):

Parameters
----------
exp : lsst.afw.image.Exposure
exp : `lsst.afw.image.Exposure`
The exposure in which to find and measure the brightest star.

doDisplay : bool
donutDiameter : `int` or `float`, optional
The expected diameter of donuts in pixels for use in the centre of
mass centroid measurement. If None is provided, the config option
is used.
doDisplay : `bool`
Display the image and found sources. A diplay object must have
been passed to the task constructor.

Returns
-------
result : lsst.pipe.base.Struct
result : `lsst.pipe.base.Struct`
Struct containing:
Whether the task ran successfully and found the object (bool)
The object's centroid (float, float)
Expand All @@ -445,18 +485,21 @@ def run(self, exp, doDisplay=False):
result.success=False.
"""
try:
result = self._run(exp=exp, doDisplay=doDisplay)
result = self._run(exp=exp, donutDiameter=donutDiameter, doDisplay=doDisplay)
return result
except Exception as e:
self.log.warn(f"Failed to find main source centroid {e}")
result = self._makeEmptyReturnStruct()
return result

def _run(self, exp, doDisplay=False):
def _run(self, exp, *, donutDiameter=None, doDisplay=False):
"""The actual run method, called by run()

Behaviour is documented in detail in the main run().
"""
if donutDiameter is None:
donutDiameter = self.config.donutDiameter

self.plateScale = exp.getWcs().getPixelScale().asArcseconds()
median = np.nanmedian(exp.image.array)
exp.image -= median # is put back later
Expand Down Expand Up @@ -510,11 +553,16 @@ def _run(self, exp, doDisplay=False):
self.checkResult(exp, brightestObjCentroid, brightestObjSrcNum,
self.config.centroidPixelPercentile)

boxSize = donutDiameter * 1.3 # allow some slack, as cutting off side of donut is very bad
centreOfMass = self._getCenterOfMass(exp, brightestObjCentroid, boxSize)

result = self._makeEmptyReturnStruct()
result.success = True
result.brightestObjCentroid = brightestObjCentroid
result.brightestObj_xXyY = (xx, yy)
result.brightestObjApFlux70 = brightestObjApFlux70
result.brightestObjApFlux25 = brightestObjApFlux25
result.medianXxYy = medianXxYy
result.brightestObjCentroidCofM = centreOfMass

return result
10 changes: 9 additions & 1 deletion tests/test_quickFrameMeasurement.py
Original file line number Diff line number Diff line change
Expand Up @@ -71,7 +71,10 @@ def testDirectCentroiding(self):
foundCentroid = result.brightestObjCentroid

dist = np.linalg.norm(np.asarray(foundCentroid) - np.asarray(trueCentroid))
distCoM = np.linalg.norm(np.asarray(result.brightestObjCentroidCofM) - np.asarray(trueCentroid))

self.assertLess(dist, self.TOLERANCE)
self.assertLess(distCoM, self.TOLERANCE)

# offset size shouldn't really matter, just make it >> PSF, and make
# sure the value isn't off-chip or right by the edge for any of the
Expand All @@ -95,9 +98,14 @@ def testDispersedCentroiding(self):
foundCentroid = result.brightestObjCentroid

dist = np.linalg.norm(np.asarray(foundCentroid) - np.asarray(trueCentroid))

self.assertLess(dist, self.TOLERANCE)

# this doesn't matter much as CoM centroid is for donuts
# and we don't intend to do dispersed donuts, but good to catch
# big regressions. Can be removed later if needed
distCoM = np.linalg.norm(np.asarray(result.brightestObjCentroidCofM) - np.asarray(trueCentroid))
self.assertLess(distCoM, self.TOLERANCE)


def setup_module(module):
lsst.utils.tests.init()
Expand Down