LSST Applications  21.0.0-147-g0e635eb1+1acddb5be5,22.0.0+052faf71bd,22.0.0+1ea9a8b2b2,22.0.0+6312710a6c,22.0.0+729191ecac,22.0.0+7589c3a021,22.0.0+9f079a9461,22.0.1-1-g7d6de66+b8044ec9de,22.0.1-1-g87000a6+536b1ee016,22.0.1-1-g8e32f31+6312710a6c,22.0.1-10-gd060f87+016f7cdc03,22.0.1-12-g9c3108e+df145f6f68,22.0.1-16-g314fa6d+c825727ab8,22.0.1-19-g93a5c75+d23f2fb6d8,22.0.1-19-gb93eaa13+aab3ef7709,22.0.1-2-g8ef0a89+b8044ec9de,22.0.1-2-g92698f7+9f079a9461,22.0.1-2-ga9b0f51+052faf71bd,22.0.1-2-gac51dbf+052faf71bd,22.0.1-2-gb66926d+6312710a6c,22.0.1-2-gcb770ba+09e3807989,22.0.1-20-g32debb5+b8044ec9de,22.0.1-23-gc2439a9a+fb0756638e,22.0.1-3-g496fd5d+09117f784f,22.0.1-3-g59f966b+1e6ba2c031,22.0.1-3-g849a1b8+f8b568069f,22.0.1-3-gaaec9c0+c5c846a8b1,22.0.1-32-g5ddfab5d3+60ce4897b0,22.0.1-4-g037fbe1+64e601228d,22.0.1-4-g8623105+b8044ec9de,22.0.1-5-g096abc9+d18c45d440,22.0.1-5-g15c806e+57f5c03693,22.0.1-7-gba73697+57f5c03693,master-g6e05de7fdc+c1283a92b8,master-g72cdda8301+729191ecac,w.2021.39
LSST Data Management Base Package
sourceDeblendTask.py
Go to the documentation of this file.
1 # This file is part of meas_deblender.
2 #
3 # Developed for the LSST Data Management System.
4 # This product includes software developed by the LSST Project
5 # (https://www.lsst.org).
6 # See the COPYRIGHT file at the top-level directory of this distribution
7 # for details of code ownership.
8 #
9 # This program is free software: you can redistribute it and/or modify
10 # it under the terms of the GNU General Public License as published by
11 # the Free Software Foundation, either version 3 of the License, or
12 # (at your option) any later version.
13 #
14 # This program is distributed in the hope that it will be useful,
15 # but WITHOUT ANY WARRANTY; without even the implied warranty of
16 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 # GNU General Public License for more details.
18 #
19 # You should have received a copy of the GNU General Public License
20 # along with this program. If not, see <https://www.gnu.org/licenses/>.
21 
22 __all__ = ['SourceDeblendConfig', 'SourceDeblendTask']
23 
24 import math
25 import numpy as np
26 
27 import lsst.log
28 import lsst.pex.config as pexConfig
29 import lsst.pipe.base as pipeBase
30 import lsst.afw.math as afwMath
31 import lsst.geom as geom
32 import lsst.afw.geom.ellipses as afwEll
33 import lsst.afw.image as afwImage
34 import lsst.afw.detection as afwDet
35 import lsst.afw.table as afwTable
36 
37 logger = lsst.log.Log.getLogger("meas.deblender.deblend")
38 
39 
40 class SourceDeblendConfig(pexConfig.Config):
41 
42  edgeHandling = pexConfig.ChoiceField(
43  doc='What to do when a peak to be deblended is close to the edge of the image',
44  dtype=str, default='ramp',
45  allowed={
46  'clip': 'Clip the template at the edge AND the mirror of the edge.',
47  'ramp': 'Ramp down flux at the image edge by the PSF',
48  'noclip': 'Ignore the edge when building the symmetric template.',
49  }
50  )
51 
52  strayFluxToPointSources = pexConfig.ChoiceField(
53  doc='When the deblender should attribute stray flux to point sources',
54  dtype=str, default='necessary',
55  allowed={
56  'necessary': 'When there is not an extended object in the footprint',
57  'always': 'Always',
58  'never': ('Never; stray flux will not be attributed to any deblended child '
59  'if the deblender thinks all peaks look like point sources'),
60  }
61  )
62 
63  assignStrayFlux = pexConfig.Field(dtype=bool, default=True,
64  doc='Assign stray flux (not claimed by any child in the deblender) '
65  'to deblend children.')
66 
67  strayFluxRule = pexConfig.ChoiceField(
68  doc='How to split flux among peaks',
69  dtype=str, default='trim',
70  allowed={
71  'r-to-peak': '~ 1/(1+R^2) to the peak',
72  'r-to-footprint': ('~ 1/(1+R^2) to the closest pixel in the footprint. '
73  'CAUTION: this can be computationally expensive on large footprints!'),
74  'nearest-footprint': ('Assign 100% to the nearest footprint (using L-1 norm aka '
75  'Manhattan distance)'),
76  'trim': ('Shrink the parent footprint to pixels that are not assigned to children')
77  }
78  )
79 
80  clipStrayFluxFraction = pexConfig.Field(dtype=float, default=0.001,
81  doc=('When splitting stray flux, clip fractions below '
82  'this value to zero.'))
83  psfChisq1 = pexConfig.Field(dtype=float, default=1.5, optional=False,
84  doc=('Chi-squared per DOF cut for deciding a source is '
85  'a PSF during deblending (un-shifted PSF model)'))
86  psfChisq2 = pexConfig.Field(dtype=float, default=1.5, optional=False,
87  doc=('Chi-squared per DOF cut for deciding a source is '
88  'PSF during deblending (shifted PSF model)'))
89  psfChisq2b = pexConfig.Field(dtype=float, default=1.5, optional=False,
90  doc=('Chi-squared per DOF cut for deciding a source is '
91  'a PSF during deblending (shifted PSF model #2)'))
92  maxNumberOfPeaks = pexConfig.Field(dtype=int, default=0,
93  doc=("Only deblend the brightest maxNumberOfPeaks peaks in the parent"
94  " (<= 0: unlimited)"))
95  maxFootprintArea = pexConfig.Field(dtype=int, default=1000000,
96  doc=("Maximum area for footprints before they are ignored as large; "
97  "non-positive means no threshold applied"))
98  maxFootprintSize = pexConfig.Field(dtype=int, default=0,
99  doc=("Maximum linear dimension for footprints before they are ignored "
100  "as large; non-positive means no threshold applied"))
101  minFootprintAxisRatio = pexConfig.Field(dtype=float, default=0.0,
102  doc=("Minimum axis ratio for footprints before they are ignored "
103  "as large; non-positive means no threshold applied"))
104  notDeblendedMask = pexConfig.Field(dtype=str, default="NOT_DEBLENDED", optional=True,
105  doc="Mask name for footprints not deblended, or None")
106 
107  tinyFootprintSize = pexConfig.RangeField(dtype=int, default=2, min=2, inclusiveMin=True,
108  doc=('Footprints smaller in width or height than this value '
109  'will be ignored; minimum of 2 due to PSF gradient '
110  'calculation.'))
111 
112  propagateAllPeaks = pexConfig.Field(dtype=bool, default=False,
113  doc=('Guarantee that all peaks produce a child source.'))
114  catchFailures = pexConfig.Field(
115  dtype=bool, default=False,
116  doc=("If True, catch exceptions thrown by the deblender, log them, "
117  "and set a flag on the parent, instead of letting them propagate up"))
118  maskPlanes = pexConfig.ListField(dtype=str, default=["SAT", "INTRP", "NO_DATA"],
119  doc="Mask planes to ignore when performing statistics")
120  maskLimits = pexConfig.DictField(
121  keytype=str,
122  itemtype=float,
123  default={},
124  doc=("Mask planes with the corresponding limit on the fraction of masked pixels. "
125  "Sources violating this limit will not be deblended."),
126  )
127  weightTemplates = pexConfig.Field(
128  dtype=bool, default=False,
129  doc=("If true, a least-squares fit of the templates will be done to the "
130  "full image. The templates will be re-weighted based on this fit."))
131  removeDegenerateTemplates = pexConfig.Field(dtype=bool, default=False,
132  doc=("Try to remove similar templates?"))
133  maxTempDotProd = pexConfig.Field(
134  dtype=float, default=0.5,
135  doc=("If the dot product between two templates is larger than this value, we consider them to be "
136  "describing the same object (i.e. they are degenerate). If one of the objects has been "
137  "labeled as a PSF it will be removed, otherwise the template with the lowest value will "
138  "be removed."))
139  medianSmoothTemplate = pexConfig.Field(dtype=bool, default=True,
140  doc="Apply a smoothing filter to all of the template images")
141 
142  # Testing options
143  # Some obs packages and ci packages run the full pipeline on a small
144  # subset of data to test that the pipeline is functioning properly.
145  # This is not meant as scientific validation, so it can be useful
146  # to only run on a small subset of the data that is large enough to
147  # test the desired pipeline features but not so long that the deblender
148  # is the tall pole in terms of execution times.
149  useCiLimits = pexConfig.Field(
150  dtype=bool, default=False,
151  doc="Limit the number of sources deblended for CI to prevent long build times")
152  ciDeblendChildRange = pexConfig.ListField(
153  dtype=int, default=[2, 10],
154  doc="Only deblend parent Footprints with a number of peaks in the (inclusive) range indicated."
155  "If `useCiLimits==False` then this parameter is ignored.")
156  ciNumParentsToDeblend = pexConfig.Field(
157  dtype=int, default=10,
158  doc="Only use the first `ciNumParentsToDeblend` parent footprints with a total peak count "
159  "within `ciDebledChildRange`. "
160  "If `useCiLimits==False` then this parameter is ignored.")
161 
162 
163 class SourceDeblendTask(pipeBase.Task):
164  """Split blended sources into individual sources.
165 
166  This task has no return value; it only modifies the SourceCatalog in-place.
167  """
168  ConfigClass = SourceDeblendConfig
169  _DefaultName = "sourceDeblend"
170 
171  def __init__(self, schema, peakSchema=None, **kwargs):
172  """Create the task, adding necessary fields to the given schema.
173 
174  Parameters
175  ----------
176  schema : `lsst.afw.table.Schema`
177  Schema object for measurement fields; will be modified in-place.
178  peakSchema : `lsst.afw.table.peakSchema`
179  Schema of Footprint Peaks that will be passed to the deblender.
180  Any fields beyond the PeakTable minimal schema will be transferred
181  to the main source Schema. If None, no fields will be transferred
182  from the Peaks
183  **kwargs
184  Additional keyword arguments passed to ~lsst.pipe.base.task
185  """
186  pipeBase.Task.__init__(self, **kwargs)
187  self.schemaschema = schema
188  self.toCopyFromParenttoCopyFromParent = [item.key for item in self.schemaschema
189  if item.field.getName().startswith("merge_footprint")]
190  peakMinimalSchema = afwDet.PeakTable.makeMinimalSchema()
191  if peakSchema is None:
192  # In this case, the peakSchemaMapper will transfer nothing, but we'll still have one
193  # to simplify downstream code
194  self.peakSchemaMapperpeakSchemaMapper = afwTable.SchemaMapper(peakMinimalSchema, schema)
195  else:
196  self.peakSchemaMapperpeakSchemaMapper = afwTable.SchemaMapper(peakSchema, schema)
197  for item in peakSchema:
198  if item.key not in peakMinimalSchema:
199  self.peakSchemaMapperpeakSchemaMapper.addMapping(item.key, item.field)
200  # Because SchemaMapper makes a copy of the output schema you give its ctor, it isn't
201  # updating this Schema in place. That's probably a design flaw, but in the meantime,
202  # we'll keep that schema in sync with the peakSchemaMapper.getOutputSchema() manually,
203  # by adding the same fields to both.
204  schema.addField(item.field)
205  assert schema == self.peakSchemaMapperpeakSchemaMapper.getOutputSchema(), "Logic bug mapping schemas"
206  self.addSchemaKeysaddSchemaKeys(schema)
207 
208  def addSchemaKeys(self, schema):
209  self.nChildKeynChildKey = schema.addField('deblend_nChild', type=np.int32,
210  doc='Number of children this object has (defaults to 0)')
211  self.psfKeypsfKey = schema.addField('deblend_deblendedAsPsf', type='Flag',
212  doc='Deblender thought this source looked like a PSF')
213  self.psfCenterKeypsfCenterKey = afwTable.Point2DKey.addFields(schema, 'deblend_psfCenter',
214  'If deblended-as-psf, the PSF centroid', "pixel")
215  self.psfFluxKeypsfFluxKey = schema.addField('deblend_psf_instFlux', type='D',
216  doc='If deblended-as-psf, the instrumental PSF flux', units='count')
217  self.tooManyPeaksKeytooManyPeaksKey = schema.addField('deblend_tooManyPeaks', type='Flag',
218  doc='Source had too many peaks; '
219  'only the brightest were included')
220  self.tooBigKeytooBigKey = schema.addField('deblend_parentTooBig', type='Flag',
221  doc='Parent footprint covered too many pixels')
222  self.maskedKeymaskedKey = schema.addField('deblend_masked', type='Flag',
223  doc='Parent footprint was predominantly masked')
224 
225  if self.config.catchFailures:
226  self.deblendFailedKeydeblendFailedKey = schema.addField('deblend_failed', type='Flag',
227  doc="Deblending failed on source")
228 
229  self.deblendSkippedKeydeblendSkippedKey = schema.addField('deblend_skipped', type='Flag',
230  doc="Deblender skipped this source")
231 
232  self.deblendRampedTemplateKeydeblendRampedTemplateKey = schema.addField(
233  'deblend_rampedTemplate', type='Flag',
234  doc=('This source was near an image edge and the deblender used '
235  '"ramp" edge-handling.'))
236 
237  self.deblendPatchedTemplateKeydeblendPatchedTemplateKey = schema.addField(
238  'deblend_patchedTemplate', type='Flag',
239  doc=('This source was near an image edge and the deblender used '
240  '"patched" edge-handling.'))
241 
242  self.hasStrayFluxKeyhasStrayFluxKey = schema.addField(
243  'deblend_hasStrayFlux', type='Flag',
244  doc=('This source was assigned some stray flux'))
245 
246  self.log.trace('Added keys to schema: %s', ", ".join(str(x) for x in (
247  self.nChildKeynChildKey, self.psfKeypsfKey, self.psfCenterKeypsfCenterKey, self.psfFluxKeypsfFluxKey,
248  self.tooManyPeaksKeytooManyPeaksKey, self.tooBigKeytooBigKey)))
249  self.peakCenterpeakCenter = afwTable.Point2IKey.addFields(schema, name="deblend_peak_center",
250  doc="Center used to apply constraints in scarlet",
251  unit="pixel")
252  self.peakIdKeypeakIdKey = schema.addField("deblend_peakId", type=np.int32,
253  doc="ID of the peak in the parent footprint. "
254  "This is not unique, but the combination of 'parent'"
255  "and 'peakId' should be for all child sources. "
256  "Top level blends with no parents have 'peakId=0'")
257  self.nPeaksKeynPeaksKey = schema.addField("deblend_nPeaks", type=np.int32,
258  doc="Number of initial peaks in the blend. "
259  "This includes peaks that may have been culled "
260  "during deblending or failed to deblend")
261  self.parentNPeaksKeyparentNPeaksKey = schema.addField("deblend_parentNPeaks", type=np.int32,
262  doc="Same as deblend_n_peaks, but the number of peaks "
263  "in the parent footprint")
264 
265  @pipeBase.timeMethod
266  def run(self, exposure, sources):
267  """Get the PSF from the provided exposure and then run deblend.
268 
269  Parameters
270  ----------
271  exposure : `lsst.afw.image.Exposure`
272  Exposure to be processed
273  sources : `lsst.afw.table.SourceCatalog`
274  SourceCatalog containing sources detected on this exposure.
275  """
276  psf = exposure.getPsf()
277  assert sources.getSchema() == self.schemaschema
278  self.deblenddeblend(exposure, sources, psf)
279 
280  def _getPsfFwhm(self, psf, bbox):
281  # It should be easier to get a PSF's fwhm;
282  # https://dev.lsstcorp.org/trac/ticket/3030
283  return psf.computeShape().getDeterminantRadius() * 2.35
284 
285  @pipeBase.timeMethod
286  def deblend(self, exposure, srcs, psf):
287  """Deblend.
288 
289  Parameters
290  ----------
291  exposure : `lsst.afw.image.Exposure`
292  Exposure to be processed
293  srcs : `lsst.afw.table.SourceCatalog`
294  SourceCatalog containing sources detected on this exposure
295  psf : `lsst.afw.detection.Psf`
296  Point source function
297 
298  Returns
299  -------
300  None
301  """
302  # Cull footprints if required by ci
303  if self.config.useCiLimits:
304  self.log.info(f"Using CI catalog limits, "
305  f"the original number of sources to deblend was {len(srcs)}.")
306  # Select parents with a number of children in the range
307  # config.ciDeblendChildRange
308  minChildren, maxChildren = self.config.ciDeblendChildRange
309  nPeaks = np.array([len(src.getFootprint().peaks) for src in srcs])
310  childrenInRange = np.where((nPeaks >= minChildren) & (nPeaks <= maxChildren))[0]
311  if len(childrenInRange) < self.config.ciNumParentsToDeblend:
312  raise ValueError("Fewer than ciNumParentsToDeblend children were contained in the range "
313  "indicated by ciDeblendChildRange. Adjust this range to include more "
314  "parents.")
315  # Keep all of the isolated parents and the first
316  # `ciNumParentsToDeblend` children
317  parents = nPeaks == 1
318  children = np.zeros((len(srcs),), dtype=bool)
319  children[childrenInRange[:self.config.ciNumParentsToDeblend]] = True
320  srcs = srcs[parents | children]
321  # We need to update the IdFactory, otherwise the the source ids
322  # will not be sequential
323  idFactory = srcs.getIdFactory()
324  maxId = np.max(srcs["id"])
325  idFactory.notify(maxId)
326 
327  self.log.info("Deblending %d sources" % len(srcs))
328 
329  from lsst.meas.deblender.baseline import deblend
330 
331  # find the median stdev in the image...
332  mi = exposure.getMaskedImage()
333  statsCtrl = afwMath.StatisticsControl()
334  statsCtrl.setAndMask(mi.getMask().getPlaneBitMask(self.config.maskPlanes))
335  stats = afwMath.makeStatistics(mi.getVariance(), mi.getMask(), afwMath.MEDIAN, statsCtrl)
336  sigma1 = math.sqrt(stats.getValue(afwMath.MEDIAN))
337  self.log.trace('sigma1: %g', sigma1)
338 
339  n0 = len(srcs)
340  nparents = 0
341  for i, src in enumerate(srcs):
342  # t0 = time.clock()
343 
344  fp = src.getFootprint()
345  pks = fp.getPeaks()
346 
347  # Since we use the first peak for the parent object, we should propagate its flags
348  # to the parent source.
349  src.assign(pks[0], self.peakSchemaMapperpeakSchemaMapper)
350 
351  if len(pks) < 2:
352  continue
353 
354  if self.isLargeFootprintisLargeFootprint(fp):
355  src.set(self.tooBigKeytooBigKey, True)
356  self.skipParentskipParent(src, mi.getMask())
357  self.log.warn('Parent %i: skipping large footprint (area: %i)',
358  int(src.getId()), int(fp.getArea()))
359  continue
360  if self.isMaskedisMasked(fp, exposure.getMaskedImage().getMask()):
361  src.set(self.maskedKeymaskedKey, True)
362  self.skipParentskipParent(src, mi.getMask())
363  self.log.warn('Parent %i: skipping masked footprint (area: %i)',
364  int(src.getId()), int(fp.getArea()))
365  continue
366 
367  nparents += 1
368  bb = fp.getBBox()
369  psf_fwhm = self._getPsfFwhm_getPsfFwhm(psf, bb)
370 
371  self.log.trace('Parent %i: deblending %i peaks', int(src.getId()), len(pks))
372 
373  self.preSingleDeblendHookpreSingleDeblendHook(exposure, srcs, i, fp, psf, psf_fwhm, sigma1)
374  npre = len(srcs)
375 
376  # This should really be set in deblend, but deblend doesn't have access to the src
377  src.set(self.tooManyPeaksKeytooManyPeaksKey, len(fp.getPeaks()) > self.config.maxNumberOfPeaks)
378 
379  try:
380  res = deblend(
381  fp, mi, psf, psf_fwhm, sigma1=sigma1,
382  psfChisqCut1=self.config.psfChisq1,
383  psfChisqCut2=self.config.psfChisq2,
384  psfChisqCut2b=self.config.psfChisq2b,
385  maxNumberOfPeaks=self.config.maxNumberOfPeaks,
386  strayFluxToPointSources=self.config.strayFluxToPointSources,
387  assignStrayFlux=self.config.assignStrayFlux,
388  strayFluxAssignment=self.config.strayFluxRule,
389  rampFluxAtEdge=(self.config.edgeHandling == 'ramp'),
390  patchEdges=(self.config.edgeHandling == 'noclip'),
391  tinyFootprintSize=self.config.tinyFootprintSize,
392  clipStrayFluxFraction=self.config.clipStrayFluxFraction,
393  weightTemplates=self.config.weightTemplates,
394  removeDegenerateTemplates=self.config.removeDegenerateTemplates,
395  maxTempDotProd=self.config.maxTempDotProd,
396  medianSmoothTemplate=self.config.medianSmoothTemplate
397  )
398  if self.config.catchFailures:
399  src.set(self.deblendFailedKeydeblendFailedKey, False)
400  except Exception as e:
401  if self.config.catchFailures:
402  self.log.warn("Unable to deblend source %d: %s" % (src.getId(), e))
403  src.set(self.deblendFailedKeydeblendFailedKey, True)
404  import traceback
405  traceback.print_exc()
406  continue
407  else:
408  raise
409 
410  kids = []
411  nchild = 0
412  for j, peak in enumerate(res.deblendedParents[0].peaks):
413  heavy = peak.getFluxPortion()
414  if heavy is None or peak.skip:
415  src.set(self.deblendSkippedKeydeblendSkippedKey, True)
416  if not self.config.propagateAllPeaks:
417  # Don't care
418  continue
419  # We need to preserve the peak: make sure we have enough info to create a minimal
420  # child src
421  self.log.trace("Peak at (%i,%i) failed. Using minimal default info for child.",
422  pks[j].getIx(), pks[j].getIy())
423  if heavy is None:
424  # copy the full footprint and strip out extra peaks
425  foot = afwDet.Footprint(src.getFootprint())
426  peakList = foot.getPeaks()
427  peakList.clear()
428  peakList.append(peak.peak)
429  zeroMimg = afwImage.MaskedImageF(foot.getBBox())
430  heavy = afwDet.makeHeavyFootprint(foot, zeroMimg)
431  if peak.deblendedAsPsf:
432  if peak.psfFitFlux is None:
433  peak.psfFitFlux = 0.0
434  if peak.psfFitCenter is None:
435  peak.psfFitCenter = (peak.peak.getIx(), peak.peak.getIy())
436 
437  assert(len(heavy.getPeaks()) == 1)
438 
439  src.set(self.deblendSkippedKeydeblendSkippedKey, False)
440  child = srcs.addNew()
441  nchild += 1
442  for key in self.toCopyFromParenttoCopyFromParent:
443  child.set(key, src.get(key))
444  child.assign(heavy.getPeaks()[0], self.peakSchemaMapperpeakSchemaMapper)
445  child.setParent(src.getId())
446  child.setFootprint(heavy)
447  child.set(self.psfKeypsfKey, peak.deblendedAsPsf)
448  child.set(self.hasStrayFluxKeyhasStrayFluxKey, peak.strayFlux is not None)
449  if peak.deblendedAsPsf:
450  (cx, cy) = peak.psfFitCenter
451  child.set(self.psfCenterKeypsfCenterKey, geom.Point2D(cx, cy))
452  child.set(self.psfFluxKeypsfFluxKey, peak.psfFitFlux)
453  child.set(self.deblendRampedTemplateKeydeblendRampedTemplateKey, peak.hasRampedTemplate)
454  child.set(self.deblendPatchedTemplateKeydeblendPatchedTemplateKey, peak.patched)
455 
456  # Set the position of the peak from the parent footprint
457  # This will make it easier to match the same source across
458  # deblenders and across observations, where the peak
459  # position is unlikely to change unless enough time passes
460  # for a source to move on the sky.
461  child.set(self.peakCenterpeakCenter, geom.Point2I(pks[j].getIx(), pks[j].getIy()))
462  child.set(self.peakIdKeypeakIdKey, pks[j].getId())
463 
464  # The children have a single peak
465  child.set(self.nPeaksKeynPeaksKey, 1)
466  # Set the number of peaks in the parent
467  child.set(self.parentNPeaksKeyparentNPeaksKey, len(pks))
468 
469  kids.append(child)
470 
471  # Child footprints may extend beyond the full extent of their parent's which
472  # results in a failure of the replace-by-noise code to reinstate these pixels
473  # to their original values. The following updates the parent footprint
474  # in-place to ensure it contains the full union of itself and all of its
475  # children's footprints.
476  spans = src.getFootprint().spans
477  for child in kids:
478  spans = spans.union(child.getFootprint().spans)
479  src.getFootprint().setSpans(spans)
480 
481  src.set(self.nChildKeynChildKey, nchild)
482 
483  self.postSingleDeblendHookpostSingleDeblendHook(exposure, srcs, i, npre, kids, fp, psf, psf_fwhm, sigma1, res)
484  # print('Deblending parent id', src.getId(), 'took', time.clock() - t0)
485 
486  n1 = len(srcs)
487  self.log.info('Deblended: of %i sources, %i were deblended, creating %i children, total %i sources'
488  % (n0, nparents, n1-n0, n1))
489 
490  def preSingleDeblendHook(self, exposure, srcs, i, fp, psf, psf_fwhm, sigma1):
491  pass
492 
493  def postSingleDeblendHook(self, exposure, srcs, i, npre, kids, fp, psf, psf_fwhm, sigma1, res):
494  pass
495 
496  def isLargeFootprint(self, footprint):
497  """Returns whether a Footprint is large
498 
499  'Large' is defined by thresholds on the area, size and axis ratio.
500  These may be disabled independently by configuring them to be non-positive.
501 
502  This is principally intended to get rid of satellite streaks, which the
503  deblender or other downstream processing can have trouble dealing with
504  (e.g., multiple large HeavyFootprints can chew up memory).
505  """
506  if self.config.maxFootprintArea > 0 and footprint.getArea() > self.config.maxFootprintArea:
507  return True
508  if self.config.maxFootprintSize > 0:
509  bbox = footprint.getBBox()
510  if max(bbox.getWidth(), bbox.getHeight()) > self.config.maxFootprintSize:
511  return True
512  if self.config.minFootprintAxisRatio > 0:
513  axes = afwEll.Axes(footprint.getShape())
514  if axes.getB() < self.config.minFootprintAxisRatio*axes.getA():
515  return True
516  return False
517 
518  def isMasked(self, footprint, mask):
519  """Returns whether the footprint violates the mask limits
520  """
521  size = float(footprint.getArea())
522  for maskName, limit in self.config.maskLimits.items():
523  maskVal = mask.getPlaneBitMask(maskName)
524  unmaskedSpan = footprint.spans.intersectNot(mask, maskVal) # spanset of unmasked pixels
525  if (size - unmaskedSpan.getArea())/size > limit:
526  return True
527  return False
528 
529  def skipParent(self, source, mask):
530  """Indicate that the parent source is not being deblended
531 
532  We set the appropriate flags and mask.
533 
534  Parameters
535  ----------
536  source : `lsst.afw.table.SourceRecord`
537  The source to flag as skipped
538  mask : `lsst.afw.image.Mask`
539  The mask to update
540  """
541  fp = source.getFootprint()
542  source.set(self.deblendSkippedKeydeblendSkippedKey, True)
543  if self.config.notDeblendedMask:
544  mask.addMaskPlane(self.config.notDeblendedMask)
545  fp.spans.setMask(mask, mask.getPlaneBitMask(self.config.notDeblendedMask))
546 
547  # Set the center of the parent
548  bbox = fp.getBBox()
549  centerX = int(bbox.getMinX()+bbox.getWidth()/2)
550  centerY = int(bbox.getMinY()+bbox.getHeight()/2)
551  source.set(self.peakCenterpeakCenter, geom.Point2I(centerX, centerY))
552  # There are no deblended children, so nChild = 0
553  source.set(self.nChildKeynChildKey, 0)
554  # But we also want to know how many peaks that we would have
555  # deblended if the parent wasn't skipped.
556  source.set(self.nPeaksKeynPeaksKey, len(fp.peaks))
557  # Top level parents are not a detected peak, so they have no peakId
558  source.set(self.peakIdKeypeakIdKey, 0)
559  # Top level parents also have no parentNPeaks
560  source.set(self.parentNPeaksKeyparentNPeaksKey, 0)
int max
Class to describe the properties of a detected object from an image.
Definition: Footprint.h:63
Pass parameters to a Statistics object.
Definition: Statistics.h:92
A mapping between the keys of two Schemas, used to copy data between them.
Definition: SchemaMapper.h:21
static Log getLogger(Log const &logger)
Definition: Log.h:772
def postSingleDeblendHook(self, exposure, srcs, i, npre, kids, fp, psf, psf_fwhm, sigma1, res)
def preSingleDeblendHook(self, exposure, srcs, i, fp, psf, psf_fwhm, sigma1)
def __init__(self, schema, peakSchema=None, **kwargs)
HeavyFootprint< ImagePixelT, MaskPixelT, VariancePixelT > makeHeavyFootprint(Footprint const &foot, lsst::afw::image::MaskedImage< ImagePixelT, MaskPixelT, VariancePixelT > const &img, HeavyFootprintCtrl const *ctrl=nullptr)
Create a HeavyFootprint with footprint defined by the given Footprint and pixel values from the given...
Backwards-compatibility support for depersisting the old Calib (FluxMag0/FluxMag0Err) objects.
Statistics makeStatistics(lsst::afw::image::Image< Pixel > const &img, lsst::afw::image::Mask< image::MaskPixel > const &msk, int const flags, StatisticsControl const &sctrl=StatisticsControl())
Handle a watered-down front-end to the constructor (no variance)
Definition: Statistics.h:359
Definition: Log.h:717