doMerge = pexConfig.Field(
dtype=bool,
default=True,
doc="Merge positive and negative diaSources with grow radius "
"set by growFootprint"
)
doForcedMeasurement = pexConfig.Field(
dtype=bool,
default=True,
doc="Force photometer diaSource locations on PVI?"
)
doAddMetrics = pexConfig.Field(
dtype=bool,
default=False,
doc="Add columns to the source table to hold analysis metrics?"
)
doSubtractBackground = pexConfig.Field(
dtype=bool,
doc="Subtract a background model from the image before detection?",
default=True,
)
doWriteBackground = pexConfig.Field(
dtype=bool,
doc="Persist the fitted background model?",
default=False,
)
subtractInitialBackground = pexConfig.ConfigurableField(
target=lsst.meas.algorithms.SubtractBackgroundTask,
doc="Task to perform intial background subtraction, before first detection pass.",
)
subtractFinalBackground = pexConfig.ConfigurableField(
target=lsst.meas.algorithms.SubtractBackgroundTask,
doc="Task to perform final background subtraction, after first detection pass.",
)
detection = pexConfig.ConfigurableField(
target=SourceDetectionTask,
doc="Final source detection for diaSource measurement",
)
streakDetection = pexConfig.ConfigurableField(
target=SourceDetectionTask,
doc="Separate source detection used only for streak masking",
)
doDeblend = pexConfig.Field(
dtype=bool,
default=False,
doc="Deblend DIASources after detection?"
)
deblend = pexConfig.ConfigurableField(
target=lsst.meas.deblender.SourceDeblendTask,
doc="Task to split blended sources into their components."
)
measurement = pexConfig.ConfigurableField(
target=DipoleFitTask,
doc="Task to measure sources on the difference image.",
)
doApCorr = lsst.pex.config.Field(
dtype=bool,
default=True,
doc="Run subtask to apply aperture corrections"
)
applyApCorr = lsst.pex.config.ConfigurableField(
target=ApplyApCorrTask,
doc="Task to apply aperture corrections"
)
forcedMeasurement = pexConfig.ConfigurableField(
target=ForcedMeasurementTask,
doc="Task to force photometer science image at diaSource locations.",
)
growFootprint = pexConfig.Field(
dtype=int,
default=2,
doc="Grow positive and negative footprints by this many pixels before merging"
)
diaSourceMatchRadius = pexConfig.Field(
dtype=float,
default=0.5,
doc="Match radius (in arcseconds) for DiaSource to Source association"
)
doSkySources = pexConfig.Field(
dtype=bool,
default=False,
doc="Generate sky sources?",
)
skySources = pexConfig.ConfigurableField(
target=SkyObjectsTask,
doc="Generate sky sources",
)
doMaskStreaks = pexConfig.Field(
dtype=bool,
default=True,
doc="Turn on streak masking",
)
maskStreaks = pexConfig.ConfigurableField(
target=MaskStreaksTask,
doc="Subtask for masking streaks. Only used if doMaskStreaks is True. "
"Adds a mask plane to an exposure, with the mask plane name set by streakMaskName.",
)
streakBinFactor = pexConfig.Field(
dtype=int,
default=4,
doc="Bin scale factor to use when rerunning detection for masking streaks. "
"Only used if doMaskStreaks is True.",
)
writeStreakInfo = pexConfig.Field(
dtype=bool,
default=False,
doc="Record the parameters of any detected streaks. For LSST, this should be turned off except for "
"development work."
)
findGlints = pexConfig.ConfigurableField(
target=FindGlintTrailsTask,
doc="Subtask for finding glint trails, usually caused by satellites or debris."
)
writeGlintInfo = pexConfig.Field(
dtype=bool,
default=True,
doc="Record the parameters of any detected glint trails."
)
setPrimaryFlags = pexConfig.ConfigurableField(
target=SetPrimaryFlagsTask,
doc="Task to add isPrimary and deblending-related flags to the catalog."
)
badSourceFlags = lsst.pex.config.ListField(
dtype=str,
doc="Sources with any of these flags set are removed before writing the output catalog.",
default=("base_PixelFlags_flag_offimage",
"base_PixelFlags_flag_interpolatedCenterAll",
"base_PixelFlags_flag_badCenterAll",
"base_PixelFlags_flag_edgeCenterAll",
"base_PixelFlags_flag_nodataCenterAll",
"base_PixelFlags_flag_saturatedCenterAll",
),
)
clearMaskPlanes = lsst.pex.config.ListField(
dtype=str,
doc="Mask planes to clear before running detection.",
default=("DETECTED", "DETECTED_NEGATIVE", "NOT_DEBLENDED", "STREAK"),
)
raiseOnBadSubtractionRatio = pexConfig.Field(
dtype=bool,
default=True,
doc="Raise an error if the ratio of power in detected footprints"
" on the difference image to the power in footprints on the science"
" image exceeds ``badSubtractionRatioThreshold``",
)
badSubtractionRatioThreshold = pexConfig.Field(
dtype=float,
default=0.2,
doc="Maximum ratio of power in footprints on the difference image to"
" the same footprints on the science image."
"Only used if ``raiseOnBadSubtractionRatio`` is set",
)
badSubtractionVariationThreshold = pexConfig.Field(
dtype=float,
default=0.4,
doc="Maximum standard deviation of the ratio of power in footprints on"
" the difference image to the same footprints on the science image."
"Only used if ``raiseOnBadSubtractionRatio`` is set",
)
raiseOnNoDiaSources = pexConfig.Field(
dtype=bool,
default=True,
doc="Raise an algorithm error if no diaSources are detected.",
)
run_sattle = pexConfig.Field(
dtype=bool,
default=False,
doc="If true, dia source bounding boxes will be sent for verification"
"to the sattle service."
)
sattle_historical = pexConfig.Field(
dtype=bool,
default=False,
doc="If re-running a pipeline that requires sattle, this should be set "
"to True. This will populate sattle's cache with the historic data "
"closest in time to the exposure."
)
idGenerator = DetectorVisitIdGeneratorConfig.make_field()
def setDefaults(self):
# Background subtraction
# Use a small binsize for the first pass to reduce detections on glints
# and extended structures. Should not affect the detectability of
# faint diaSources
self.subtractInitialBackground.binSize = 8
self.subtractInitialBackground.useApprox = False
self.subtractInitialBackground.statisticsProperty = "MEDIAN"
self.subtractInitialBackground.doFilterSuperPixels = True
self.subtractInitialBackground.ignoredPixelMask = ["BAD",
"EDGE",
"DETECTED",
"DETECTED_NEGATIVE",
"NO_DATA",
]
# Use a larger binsize for the final background subtraction, to reduce
# over-subtraction of bright objects.
self.subtractFinalBackground.binSize = 40
self.subtractFinalBackground.useApprox = False
self.subtractFinalBackground.statisticsProperty = "MEDIAN"
self.subtractFinalBackground.doFilterSuperPixels = True
self.subtractFinalBackground.ignoredPixelMask = ["BAD",
"EDGE",
"DETECTED",
"DETECTED_NEGATIVE",
"NO_DATA",
]
# DiaSource Detection
self.detection.thresholdPolarity = "both"
self.detection.thresholdValue = 5.0
self.detection.reEstimateBackground = False
self.detection.thresholdType = "pixel_stdev"
self.detection.excludeMaskPlanes = ["EDGE",
"BAD",
]
# Copy configs for binned streak detection from the base detection task
self.streakDetection.thresholdType = self.detection.thresholdType
self.streakDetection.reEstimateBackground = False
self.streakDetection.excludeMaskPlanes = self.detection.excludeMaskPlanes
self.streakDetection.thresholdValue = self.detection.thresholdValue
# Only detect positive streaks
self.streakDetection.thresholdPolarity = "positive"
# Do not grow detected mask for streaks
self.streakDetection.nSigmaToGrow = 0
# Set the streak mask along the entire fit line, not only where the
# detected mask is set.
self.maskStreaks.onlyMaskDetected = False
# Restrict streak masking from growing too large
self.maskStreaks.maxStreakWidth = 100
# Restrict the number of iterations allowed for fitting streaks
# When the fit is good it should solve quickly, and exit a bad fit quickly
self.maskStreaks.maxFitIter = 10
# Only mask to 2 sigma in width
self.maskStreaks.nSigmaMask = 2
# Threshold for including streaks after the Hough Transform.
# A lower value will detect more features that are less linear.
self.maskStreaks.absMinimumKernelHeight = 2
self.measurement.plugins.names |= ["ext_trailedSources_Naive",
"base_LocalPhotoCalib",
"base_LocalWcs",
"ext_shapeHSM_HsmSourceMoments",
"ext_shapeHSM_HsmPsfMoments",
"base_ClassificationSizeExtendedness",
]
self.measurement.slots.psfShape = "ext_shapeHSM_HsmPsfMoments"
self.measurement.slots.shape = "ext_shapeHSM_HsmSourceMoments"
self.measurement.plugins["base_SdssCentroid"].maxDistToPeak = 5.0
self.forcedMeasurement.plugins = ["base_TransformedCentroid", "base_PsfFlux"]
self.forcedMeasurement.copyColumns = {
"id": "objectId", "parent": "parentObjectId", "coord_ra": "coord_ra", "coord_dec": "coord_dec"}
self.forcedMeasurement.slots.centroid = "base_TransformedCentroid"
self.forcedMeasurement.slots.shape = None
# Keep track of which footprints contain streaks
self.measurement.plugins["base_PixelFlags"].masksFpAnywhere = [
"STREAK", "INJECTED", "INJECTED_TEMPLATE"]
self.measurement.plugins["base_PixelFlags"].masksFpCenter = [
"STREAK", "INJECTED", "INJECTED_TEMPLATE"]
self.skySources.avoidMask = ["DETECTED", "DETECTED_NEGATIVE", "BAD", "NO_DATA", "EDGE"]
def validate(self):
super().validate()
if self.run_sattle:
if not os.getenv("SATTLE_URI_BASE"):
raise pexConfig.FieldValidationError(DetectAndMeasureConfig.run_sattle, self,
"Sattle requested but SATTLE_URI_BASE "
"environment variable not set.")
class DetectAndMeasureTask(lsst.pipe.base.PipelineTask):
ConfigClass = DetectAndMeasureConfig
_DefaultName = "detectAndMeasure"
def __init__(self, **kwargs):
super().__init__(**kwargs)
self.schema = afwTable.SourceTable.makeMinimalSchema()
self.algMetadata = dafBase.PropertyList()
if self.config.doSubtractBackground:
self.makeSubtask("subtractInitialBackground")
self.makeSubtask("subtractFinalBackground")
self.makeSubtask("detection", schema=self.schema)
if self.config.doDeblend:
self.makeSubtask("deblend", schema=self.schema)
self.makeSubtask("setPrimaryFlags", schema=self.schema, isSingleFrame=True)
self.makeSubtask("measurement", schema=self.schema,
algMetadata=self.algMetadata)
if self.config.doApCorr:
self.makeSubtask("applyApCorr", schema=self.measurement.schema)
if self.config.doForcedMeasurement:
self.schema.addField(
"ip_diffim_forced_PsfFlux_instFlux", "D",
"Forced PSF flux measured on the direct image.",
units="count")
self.schema.addField(
"ip_diffim_forced_PsfFlux_instFluxErr", "D",
"Forced PSF flux error measured on the direct image.",
units="count")
self.schema.addField(
"ip_diffim_forced_PsfFlux_area", "F",
"Forced PSF flux effective area of PSF.",
units="pixel")
self.schema.addField(
"ip_diffim_forced_PsfFlux_flag", "Flag",
"Forced PSF flux general failure flag.")
self.schema.addField(
"ip_diffim_forced_PsfFlux_flag_noGoodPixels", "Flag",
"Forced PSF flux not enough non-rejected pixels in data to attempt the fit.")
self.schema.addField(
"ip_diffim_forced_PsfFlux_flag_edge", "Flag",
"Forced PSF flux object was too close to the edge of the image to use the full PSF model.")
self.schema.addField(
"ip_diffim_forced_template_PsfFlux_instFlux", "D",
"Forced PSF flux measured on the template image.",
units="count")
self.schema.addField(
"ip_diffim_forced_template_PsfFlux_instFluxErr", "D",
"Forced PSF flux error measured on the template image.",
units="count")
self.schema.addField(
"ip_diffim_forced_template_PsfFlux_area", "F",
"Forced template PSF flux effective area of PSF.",
units="pixel")
self.schema.addField(
"ip_diffim_forced_template_PsfFlux_flag", "Flag",
"Forced template PSF flux general failure flag.")
self.schema.addField(
"ip_diffim_forced_template_PsfFlux_flag_noGoodPixels", "Flag",
"Forced template PSF flux not enough non-rejected pixels in data to attempt the fit.")
self.schema.addField(
"ip_diffim_forced_template_PsfFlux_flag_edge", "Flag",
)
self.makeSubtask("forcedMeasurement", refSchema=self.schema)
self.schema.addField("refMatchId", "L", "unique id of reference catalog match")
self.schema.addField("srcMatchId", "L", "unique id of source match")
# Create the sky source task for use by metrics,
# even if sky sources are not added to the diaSource catalog
self.makeSubtask("skySources", schema=self.schema)
if self.config.doMaskStreaks:
self.makeSubtask("maskStreaks")
self.makeSubtask("streakDetection")
self.makeSubtask("findGlints")
self.schema.addField("glint_trail", "Flag", "DiaSource is part of a glint trail.")
self.schema.addField("reliability", type="F", doc="Reliability score of the DiaSource")
# To get the "merge_*" fields in the schema; have to re-initialize
# this later, once we have a peak schema post-detection.
lsst.afw.detection.FootprintMergeList(self.schema, ["positive", "negative"])
# Check that the schema and config are consistent
for flag in self.config.badSourceFlags:
if flag not in self.schema:
raise pipeBase.InvalidQuantumError("Field %s not in schema" % flag)
# initialize InitOutputs
self.outputSchema = afwTable.SourceCatalog(self.schema)
self.outputSchema.getTable().setMetadata(self.algMetadata)
def runQuantum(self, butlerQC: pipeBase.QuantumContext,
inputRefs: pipeBase.InputQuantizedConnection,
outputRefs: pipeBase.OutputQuantizedConnection):
inputs = butlerQC.get(inputRefs)
idGenerator = self.config.idGenerator.apply(butlerQC.quantum.dataId)
idFactory = idGenerator.make_table_id_factory()
# Specify the fields that `annotate` needs below, to ensure they
# exist, even as None.
measurementResults = pipeBase.Struct(
subtractedMeasuredExposure=None,
diaSources=None,
maskedStreaks=None,
differenceBackground=None,
)
try:
self.run(**inputs, idFactory=idFactory, measurementResults=measurementResults)
except pipeBase.AlgorithmError as e:
error = pipeBase.AnnotatedPartialOutputsError.annotate(
e,
self,
measurementResults.subtractedMeasuredExposure,
measurementResults.diaSources,
measurementResults.maskedStreaks,
log=self.log
)
butlerQC.put(measurementResults, outputRefs)
raise error from e
butlerQC.put(measurementResults, outputRefs)
@timeMethod
def run(self, science, matchedTemplate, difference, kernelSources,
idFactory=None, measurementResults=None):
# Check that we have a valid PSF now before we do more work
sigma = difference.psf.computeShape(difference.psf.getAveragePosition()).getDeterminantRadius()
if np.isnan(sigma):
raise pipeBase.UpstreamFailureNoWorkFound("Invalid PSF detected! PSF width evaluates to NaN.")
# Ensure that we start with an empty detection and deblended mask.
mask = difference.mask
for mp in self.config.clearMaskPlanes:
if mp not in mask.getMaskPlaneDict():
mask.addMaskPlane(mp)
mask &= ~mask.getPlaneBitMask(self.config.clearMaskPlanes)
def processResults(self, science, matchedTemplate, difference, sources, idFactory, kernelSources,
positives=None, negatives=None, measurementResults=None):
# Run forced psf photometry on the image at the diaSource locations.
# Copy the measured flux and error into the diaSource.
forcedSources = self.forcedMeasurement.generateMeasCat(image, diaSources, wcs)
self.forcedMeasurement.run(forcedSources, image, diaSources, wcs)
if template:
base_key = 'ip_diffim_forced_template_PsfFlux'
else:
base_key = 'ip_diffim_forced_PsfFlux'
mapper = afwTable.SchemaMapper(forcedSources.schema, diaSources.schema)
mapper.addMapping(forcedSources.schema.find("base_PsfFlux_instFlux")[0],
f"{base_key}_instFlux", True)
mapper.addMapping(forcedSources.schema.find("base_PsfFlux_instFluxErr")[0],
f"{base_key}_instFluxErr", True)
mapper.addMapping(forcedSources.schema.find("base_PsfFlux_area")[0],
f"{base_key}_area", True)
mapper.addMapping(forcedSources.schema.find("base_PsfFlux_flag")[0],
f"{base_key}_flag", True)
mapper.addMapping(forcedSources.schema.find("base_PsfFlux_flag_noGoodPixels")[0],
f"{base_key}_flag_noGoodPixels", True)
mapper.addMapping(forcedSources.schema.find("base_PsfFlux_flag_edge")[0],
f"{base_key}_flag_edge", True)
for diaSource, forcedSource in zip(diaSources, forcedSources):
diaSource.assign(forcedSource, mapper)
def calculateMetrics(self, science, difference, diaSources, kernelSources):
maskedImage = difference.maskedImage
# Bin the diffim to enhance low surface brightness streaks
binnedMaskedImage = afwMath.binImage(maskedImage,
self.config.streakBinFactor,
self.config.streakBinFactor)
binnedExposure = afwImage.ExposureF(binnedMaskedImage.getBBox())
binnedExposure.setMaskedImage(binnedMaskedImage)
# Clear the DETECTED mask plane before streak detection
binnedExposure.mask &= ~binnedExposure.mask.getPlaneBitMask('DETECTED')
# Rerun detection to set the DETECTED mask plane on binnedExposure
sigma = difference.psf.computeShape(difference.psf.getAveragePosition()).getDeterminantRadius()
_table = afwTable.SourceTable.make(afwTable.SourceTable.makeMinimalSchema())
self.streakDetection.run(table=_table, exposure=binnedExposure, doSmooth=True,
sigma=sigma/self.config.streakBinFactor)
binnedDetectedMaskPlane = binnedExposure.mask.array & binnedExposure.mask.getPlaneBitMask('DETECTED')
rescaledDetectedMaskPlane = binnedDetectedMaskPlane.repeat(self.config.streakBinFactor,
axis=0).repeat(self.config.streakBinFactor,
axis=1)
# Create new version of a diffim with DETECTED based on binnedExposure
streakMaskedImage = maskedImage.clone()
ysize, xsize = rescaledDetectedMaskPlane.shape
streakMaskedImage.mask.array[:ysize, :xsize] |= rescaledDetectedMaskPlane
# Detect streaks on this new version of the diffim
streaks = self.maskStreaks.run(streakMaskedImage)
streakMaskPlane = streakMaskedImage.mask.array & streakMaskedImage.mask.getPlaneBitMask('STREAK')
# Apply the new STREAK mask to the original diffim
maskedImage.mask.array |= streakMaskPlane
if self.config.writeStreakInfo:
rhos = np.array([line.rho for line in streaks.lines])
thetas = np.array([line.theta for line in streaks.lines])
sigmas = np.array([line.sigma for line in streaks.lines])
chi2s = np.array([line.reducedChi2 for line in streaks.lines])
modelMaximums = np.array([line.modelMaximum for line in streaks.lines])
streakInfo = {'rho': rhos, 'theta': thetas, 'sigma': sigmas, 'reducedChi2': chi2s,
'modelMaximum': modelMaximums}
else:
streakInfo = {'rho': np.array([]), 'theta': np.array([]), 'sigma': np.array([]),
'reducedChi2': np.array([]), 'modelMaximum': np.array([])}
return pipeBase.Struct(maskedStreaks=streakInfo)
class DetectAndMeasureScoreConnections(DetectAndMeasureConnections):
scoreExposure = pipeBase.connectionTypes.Input(
doc="Maximum likelihood image for detection.",
dimensions=("instrument", "visit", "detector"),
storageClass="ExposureF",
name="{fakesType}{coaddName}Diff_scoreExp",
)
class DetectAndMeasureScoreConfig(DetectAndMeasureConfig,
pipelineConnections=DetectAndMeasureScoreConnections):
pass
class DetectAndMeasureScoreTask(DetectAndMeasureTask):
ConfigClass = DetectAndMeasureScoreConfig
_DefaultName = "detectAndMeasureScore"
@timeMethod
def run(self, science, matchedTemplate, difference, scoreExposure, kernelSources,
idFactory=None):
Definition at line 571 of file detectAndMeasure.py.