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.

Sign up for GitHub

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-46141: Use ScienceSourceSelector to select kernel candidates #340

Merged
merged 5 commits into from
Sep 18, 2024
Merged

DM-46141: Use ScienceSourceSelector to select kernel candidates #340

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
0611d7b
Add coordinate errors to test schemas
isullivan Sep 10, 2024
ac59a76
Refactor subtraction unit tests to configure Task in common method
isullivan Sep 10, 2024
ff1e364
Switch to ScienceSourceSelectorTask
mrawls Sep 9, 2024
eed2622
Add constraints for isolated stars to psf match source selector
isullivan Sep 10, 2024
0b9d46f
Extend template mask check to neighboring pixels
isullivan Sep 10, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
64 changes: 41 additions & 23 deletions python/lsst/ip/diffim/subtractImages.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,7 @@
import lsst.afw.math
import lsst.geom
from lsst.ip.diffim.utils import evaluateMeanPsfFwhm, getPsfFwhm
from lsst.meas.algorithms import ScaleVarianceTask
from lsst.meas.algorithms import ScaleVarianceTask, ScienceSourceSelectorTask
import lsst.pex.config
import lsst.pipe.base
import lsst.pex.exceptions
Expand Down Expand Up @@ -173,6 +173,10 @@ class AlardLuptonSubtractBaseConfig(lsst.pex.config.Config):
dtype=bool,
default=False,
)
sourceSelector = lsst.pex.config.ConfigurableField(
target=ScienceSourceSelectorTask,
doc="Task to select sources to be used for PSF matching.",
)
detectionThreshold = lsst.pex.config.Field(
dtype=float,
default=10,
Expand Down Expand Up @@ -206,6 +210,8 @@ class AlardLuptonSubtractBaseConfig(lsst.pex.config.Config):
"base_PixelFlags_flag_saturated",
"base_PixelFlags_flag_bad",
),
deprecated="This field is no longer used and will be removed after v28."
"Set the equivalent field for the sourceSelector subtask instead.",
)
excludeMaskPlanes = lsst.pex.config.ListField(
dtype=str,
Expand Down Expand Up @@ -240,6 +246,13 @@ def setDefaults(self):
self.makeKernel.kernel.active.fitForBackground = self.doSubtractBackground
self.makeKernel.kernel.active.spatialKernelOrder = 1
self.makeKernel.kernel.active.spatialBgOrder = 2
self.sourceSelector.doUnresolved = True # apply star-galaxy separation
self.sourceSelector.doIsolated = True # apply isolated star selection
self.sourceSelector.doRequirePrimary = True # apply primary flag selection
self.sourceSelector.doSkySources = False # Do not include sky sources
self.sourceSelector.doSignalToNoise = True # apply signal to noise filter
self.sourceSelector.signalToNoise.minimum = 10
self.sourceSelector.signalToNoise.maximum = 500


class AlardLuptonSubtractConfig(AlardLuptonSubtractBaseConfig, lsst.pipe.base.PipelineTaskConfig,
Expand All @@ -265,6 +278,7 @@ def __init__(self, **kwargs):
super().__init__(**kwargs)
self.makeSubtask("decorrelate")
self.makeSubtask("makeKernel")
self.makeSubtask("sourceSelector")
if self.config.doScaleVariance:
self.makeSubtask("scaleVariance")

Expand Down Expand Up @@ -772,26 +786,23 @@ def _sourceSelector(self, sources, mask):
If there are too few sources to compute the PSF matching kernel
remaining after source selection.
"""
flags = np.ones(len(sources), dtype=bool)
for flag in self.config.badSourceFlags:
try:
flags *= ~sources[flag]
except Exception as e:
self.log.warning("Could not apply source flag: %s", e)
signalToNoise = sources.getPsfInstFlux()/sources.getPsfInstFluxErr()
sToNFlag = signalToNoise > self.config.detectionThreshold
flags *= sToNFlag
sToNFlagMax = signalToNoise < self.config.detectionThresholdMax
flags *= sToNFlagMax
flags *= self._checkMask(mask, sources, self.config.excludeMaskPlanes)
selectSources = sources[flags].copy(deep=True)

selected = self.sourceSelector.selectSources(sources).selected
nInitialSelected = np.count_nonzero(selected)
selected *= self._checkMask(mask, sources, self.config.excludeMaskPlanes)
nSelected = np.count_nonzero(selected)
self.log.info("Rejecting %i candidate sources: an excluded template mask plane is set.",
nInitialSelected - nSelected)
selectSources = sources[selected].copy(deep=True)
# Trim selectSources if they exceed ``maxKernelSources``.
# Keep the highest signal-to-noise sources of those selected.
if (len(selectSources) > self.config.maxKernelSources) & (self.config.maxKernelSources > 0):
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Consider adding a one-line comment along the lines of # Trim selectSources so it doesn't exceed maxKernalSources. I presume selectSources is not ordered in any meaningful way? So we won't, for example, throw out all the sources in one region of the image, or those with certain S/N or brightness properties?

signalToNoise = selectSources.getPsfInstFlux()/selectSources.getPsfInstFluxErr()
indices = np.argsort(signalToNoise)
indices = indices[-self.config.maxKernelSources:]
flags = np.zeros(len(selectSources), dtype=bool)
flags[indices] = True
selectSources = selectSources[flags].copy(deep=True)
selected = np.zeros(len(selectSources), dtype=bool)
selected[indices] = True
selectSources = selectSources[selected].copy(deep=True)

self.log.info("%i/%i=%.1f%% of sources selected for PSF matching from the input catalog",
len(selectSources), len(sources), 100*len(selectSources)/len(sources))
Expand All @@ -805,8 +816,8 @@ def _sourceSelector(self, sources, mask):
return selectSources

@staticmethod
def _checkMask(mask, sources, excludeMaskPlanes):
"""Exclude sources that are located on masked pixels.
def _checkMask(mask, sources, excludeMaskPlanes, checkAdjacent=True):
"""Exclude sources that are located on or adjacent to masked pixels.

Parameters
----------
Expand All @@ -831,11 +842,18 @@ def _checkMask(mask, sources, excludeMaskPlanes):

excludePixelMask = mask.getPlaneBitMask(setExcludeMaskPlanes)

xv = np.rint(sources.getX() - mask.getX0())
yv = np.rint(sources.getY() - mask.getY0())
xv = (np.rint(sources.getX() - mask.getX0())).astype(int)
yv = (np.rint(sources.getY() - mask.getY0())).astype(int)

mv = mask.array[yv.astype(int), xv.astype(int)]
flags = np.bitwise_and(mv, excludePixelMask) == 0
flags = np.ones(len(sources), dtype=bool)
if checkAdjacent:
pixRange = (0, -1, 1)
else:
pixRange = (0,)
for j in pixRange:
for i in pixRange:
mv = mask.array[yv + j, xv + i]
flags *= np.bitwise_and(mv, excludePixelMask) == 0
return flags

def _prepareInputs(self, template, science, visitSummary=None):
Expand Down
14 changes: 14 additions & 0 deletions python/lsst/ip/diffim/utils.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1095,6 +1095,8 @@ def detectTestSources(exposure, addMaskPlanes=None):
addMaskPlanes = ["STREAK", "INJECTED", "INJECTED_TEMPLATE"]

schema = afwTable.SourceTable.makeMinimalSchema()
# Add coordinate error fields:
afwTable.CoordKey.addErrorFields(schema)
selectDetection = measAlg.SourceDetectionTask(schema=schema)
selectMeasurement = measBase.SingleFrameMeasurementTask(schema=schema)
table = afwTable.SourceTable.make(schema)
Expand Down Expand Up @@ -1267,6 +1269,8 @@ def _makeTruthSchema():
Calib, Ap, and Psf flux slots all are set to ``truth_instFlux``.
"""
schema = afwTable.SourceTable.makeMinimalSchema()
# Add coordinate error fields:
afwTable.CoordKey.addErrorFields(schema)
keys = {}
# Don't use a FluxResultKey so we can manage the flux and err separately.
keys["instFlux"] = schema.addField("truth_instFlux", type=np.float64,
Expand All @@ -1280,6 +1284,11 @@ def _makeTruthSchema():
schema.addField("base_PixelFlags_flag_interpolated", "Flag", "testing flag.")
schema.addField("base_PixelFlags_flag_saturated", "Flag", "testing flag.")
schema.addField("base_PixelFlags_flag_bad", "Flag", "testing flag.")
schema.addField("base_PixelFlags_flag_edge", "Flag", "testing flag.")
schema.addField("base_PsfFlux_flag", "Flag", "testing flag.")
schema.addField("base_ClassificationSizeExtendedness_value", "Flag", "testing flag.")
schema.addField("deblend_nChild", "Flag", "testing flag.")
schema.addField("detect_isPrimary", "Flag", "testing flag.")
schema.getAliasMap().set("slot_Centroid", "truth")
schema.getAliasMap().set("slot_CalibFlux", "truth")
schema.getAliasMap().set("slot_ApFlux", "truth")
Expand Down Expand Up @@ -1316,6 +1325,11 @@ def _fillTruthCatalog(injectList):
footprint.addPeak(x, y, flux)
record.setFootprint(footprint)

# Set source records for isolated stars
record["base_ClassificationSizeExtendedness_value"] = 0
record["deblend_nChild"] = 0
record["detect_isPrimary"] = True

return catalog


Expand Down
2 changes: 2 additions & 0 deletions tests/test_detectAndMeasure.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -488,6 +488,8 @@ def test_fake_mask_plane_propagation(self):
template_fake_masked = (template.mask.array & template_fake_bitmask) > 0

subtractConfig = subtractImages.AlardLuptonSubtractTask.ConfigClass()
subtractConfig.sourceSelector.signalToNoise.fluxField = "truth_instFlux"
subtractConfig.sourceSelector.signalToNoise.errField = "truth_instFluxErr"
subtractTask = subtractImages.AlardLuptonSubtractTask(config=subtractConfig)
subtraction = subtractTask.run(template, science, sources)

Expand Down
Loading
Loading