LSST Applications 29.1.0,g0fba68d861+6b120c4394,g123d84c11c+8c5ae1fdc5,g1ec0fe41b4+191117f6ec,g1fd858c14a+c8450ae71a,g3533f9d6cb+a04f9ee0ab,g35bb328faa+8c5ae1fdc5,g3f0dcc2b1b+7df08700bd,g4178042926+b4254969db,g44ba364a48+04455b336b,g53246c7159+8c5ae1fdc5,g60b5630c4e+a04f9ee0ab,g663da51e9b+b05e6e1875,g67b6fd64d1+250bf6acd3,g78460c75b0+7e33a9eb6d,g786e29fd12+668abc6043,g8352419a5c+8c5ae1fdc5,g87e3079a85+d3fa38de54,g8852436030+cd899e2626,g89139ef638+250bf6acd3,g93a033419f+31ead11197,g989de1cb63+250bf6acd3,g9f33ca652e+f6053ecf14,ga1e959baac+5fbc491aed,ga2f891cd6c+a04f9ee0ab,gabe3b4be73+8856018cbb,gabf8522325+1f7e6d67b9,gac2eed3f23+250bf6acd3,gb1101e3267+0c331e9486,gb89ab40317+250bf6acd3,gcf25f946ba+cd899e2626,gd107969129+8964d67276,gd6cbbdb0b4+6bbecc8878,gde0f65d7ad+d65f9e019a,ge278dab8ac+eb3bbeb12f,ge410e46f29+250bf6acd3,gf5e32f922b+8c5ae1fdc5,gff02db199a+747430a128,gffe7e49bb4+a04f9ee0ab
LSST Data Management Base Package
Loading...
Searching...
No Matches
lsst.ip.diffim.getTemplate Namespace Reference

Classes

class  GetDcrTemplateConnections
 
class  GetTemplateConnections
 

Functions

 run (self, *, coaddExposureHandles, bbox, wcs, dataIds, physical_filter)
 
 _makeExposureCatalog (self, exposureRefs, dataIds)
 
 _merge (self, maskedImages, bbox, wcs)
 
 _makePsf (self, template, catalog, wcs)
 
 checkPatchList (self, patchList)
 
 getDcrModel (self, patchList, coaddRefs, visitInfo)
 
 _selectDataRef (coaddRef, tract, patch)
 

Variables

 detectorPolygon = geom.Box2D(bbox)
 
int overlappingArea = 0
 
 coaddExposures = collections.defaultdict(list)
 
 dataIds = collections.defaultdict(list)
 
 dataId = coaddRef.dataId
 
 patchWcs = skymap[dataId["tract"]].getWcs()
 
 patchBBox = skymap[dataId["tract"]][dataId["patch"]].getOuterBBox()
 
 patchCorners = patchWcs.pixelToSky(geom.Box2D(patchBBox).getCorners())
 
 patchPolygon = afwGeom.Polygon(wcs.skyToPixel(patchCorners))
 
 patchList = dict()
 
 subfilter = dataId["subfilter"]
 

Function Documentation

◆ _makeExposureCatalog()

lsst.ip.diffim.getTemplate._makeExposureCatalog ( self,
exposureRefs,
dataIds )
protected
Make an exposure catalog for one tract.

Parameters
----------
exposureRefs : `list` of [`lsst.daf.butler.DeferredDatasetHandle` of \
                `lsst.afw.image.Exposure`]
    Exposures to include in the catalog.
dataIds : `list` [`lsst.daf.butler.DataCoordinate`]
    Data ids of each of the included exposures; must have "tract" and
    "patch" entries.

Returns
-------
images : `dict` [`lsst.afw.image.MaskedImage`]
    MaskedImages of each of the input exposures, for warping.
catalog : `lsst.afw.table.ExposureCatalog`
    Catalog of metadata for each exposure
totalBox : `lsst.geom.Box2I`
    The union of the bounding boxes of all the input exposures.

Definition at line 408 of file getTemplate.py.

408 def _makeExposureCatalog(self, exposureRefs, dataIds):
409 """Make an exposure catalog for one tract.
410
411 Parameters
412 ----------
413 exposureRefs : `list` of [`lsst.daf.butler.DeferredDatasetHandle` of \
414 `lsst.afw.image.Exposure`]
415 Exposures to include in the catalog.
416 dataIds : `list` [`lsst.daf.butler.DataCoordinate`]
417 Data ids of each of the included exposures; must have "tract" and
418 "patch" entries.
419
420 Returns
421 -------
422 images : `dict` [`lsst.afw.image.MaskedImage`]
423 MaskedImages of each of the input exposures, for warping.
424 catalog : `lsst.afw.table.ExposureCatalog`
425 Catalog of metadata for each exposure
426 totalBox : `lsst.geom.Box2I`
427 The union of the bounding boxes of all the input exposures.
428 """
429 catalog = afwTable.ExposureCatalog(self.schema)
430 catalog.reserve(len(exposureRefs))
431 exposures = (exposureRef.get() for exposureRef in exposureRefs)
432 images = {}
433 totalBox = geom.Box2I()
434
435 for coadd, dataId in zip(exposures, dataIds):
436 images[dataId] = coadd.maskedImage
437 bbox = coadd.getBBox()
438 totalBox = totalBox.expandedTo(bbox)
439 record = catalog.addNew()
440 record.setPsf(coadd.psf)
441 record.setWcs(coadd.wcs)
442 record.setPhotoCalib(coadd.photoCalib)
443 record.setBBox(bbox)
444 record.setValidPolygon(afwGeom.Polygon(geom.Box2D(bbox).getCorners()))
445 record.set("tract", dataId["tract"])
446 record.set("patch", dataId["patch"])
447 # Weight is used by CoaddPsf, but the PSFs from overlapping patches
448 # should be very similar, so this value mostly shouldn't matter.
449 record.set("weight", 1)
450
451 return images, catalog, totalBox
452
Cartesian polygons.
Definition Polygon.h:59
A floating-point coordinate rectangle geometry.
Definition Box.h:413
An integer coordinate rectangle.
Definition Box.h:55

◆ _makePsf()

lsst.ip.diffim.getTemplate._makePsf ( self,
template,
catalog,
wcs )
protected
Return a PSF containing the PSF at each of the input regions.

Note that although this includes all the exposures from the catalog,
the PSF knows which part of the template the inputs came from, so when
evaluated at a given position it will not include inputs that never
went in to those pixels.

Parameters
----------
template : `lsst.afw.image.Exposure`
    Generated template the PSF is for.
catalog : `lsst.afw.table.ExposureCatalog`
    Catalog of exposures that went into the template that contains all
    of the input PSFs.
wcs : `lsst.afw.geom.SkyWcs`
    WCS of the template, to warp the PSFs to.

Returns
-------
coaddPsf : `lsst.meas.algorithms.CoaddPsf`
    The meta-psf constructed from all of the input catalogs.

Definition at line 527 of file getTemplate.py.

527 def _makePsf(self, template, catalog, wcs):
528 """Return a PSF containing the PSF at each of the input regions.
529
530 Note that although this includes all the exposures from the catalog,
531 the PSF knows which part of the template the inputs came from, so when
532 evaluated at a given position it will not include inputs that never
533 went in to those pixels.
534
535 Parameters
536 ----------
537 template : `lsst.afw.image.Exposure`
538 Generated template the PSF is for.
539 catalog : `lsst.afw.table.ExposureCatalog`
540 Catalog of exposures that went into the template that contains all
541 of the input PSFs.
542 wcs : `lsst.afw.geom.SkyWcs`
543 WCS of the template, to warp the PSFs to.
544
545 Returns
546 -------
547 coaddPsf : `lsst.meas.algorithms.CoaddPsf`
548 The meta-psf constructed from all of the input catalogs.
549 """
550 # CoaddPsf centroid not only must overlap image, but must overlap the
551 # part of image with data. Use centroid of region with data.
552 boolmask = template.mask.array & template.mask.getPlaneBitMask("NO_DATA") == 0
553 maskx = afwImage.makeMaskFromArray(boolmask.astype(afwImage.MaskPixel))
554 centerCoord = afwGeom.SpanSet.fromMask(maskx, 1).computeCentroid()
555
556 ctrl = self.config.coaddPsf.makeControl()
557 coaddPsf = CoaddPsf(
558 catalog, wcs, centerCoord, ctrl.warpingKernelName, ctrl.cacheSize
559 )
560 return coaddPsf
561
562

◆ _merge()

lsst.ip.diffim.getTemplate._merge ( self,
maskedImages,
bbox,
wcs )
protected
Merge the images that came from one tract into one larger image,
ignoring NaN pixels and non-finite variance pixels from individual
exposures.

Parameters
----------
maskedImages : `dict` [`lsst.afw.image.MaskedImage` or
                       `lsst.afw.image.Exposure`]
    Images to be merged into one larger bounding box.
bbox : `lsst.geom.Box2I`
    Bounding box defining the image to merge into.
wcs : `lsst.afw.geom.SkyWcs`
    WCS of all of the input images to set on the output image.

Returns
-------
merged : `lsst.afw.image.MaskedImage`
    Merged image with all of the inputs at their respective bbox
    positions.
count : `int`
    Count of the number of good pixels (those with positive weights)
    in the merged image.
included : `list` [`int`]
    List of indexes of patches that were included in the merged
    result, to be used to trim the exposure catalog.

Definition at line 453 of file getTemplate.py.

453 def _merge(self, maskedImages, bbox, wcs):
454 """Merge the images that came from one tract into one larger image,
455 ignoring NaN pixels and non-finite variance pixels from individual
456 exposures.
457
458 Parameters
459 ----------
460 maskedImages : `dict` [`lsst.afw.image.MaskedImage` or
461 `lsst.afw.image.Exposure`]
462 Images to be merged into one larger bounding box.
463 bbox : `lsst.geom.Box2I`
464 Bounding box defining the image to merge into.
465 wcs : `lsst.afw.geom.SkyWcs`
466 WCS of all of the input images to set on the output image.
467
468 Returns
469 -------
470 merged : `lsst.afw.image.MaskedImage`
471 Merged image with all of the inputs at their respective bbox
472 positions.
473 count : `int`
474 Count of the number of good pixels (those with positive weights)
475 in the merged image.
476 included : `list` [`int`]
477 List of indexes of patches that were included in the merged
478 result, to be used to trim the exposure catalog.
479 """
480 merged = afwImage.ExposureF(bbox, wcs)
481 weights = afwImage.ImageF(bbox)
482 included = [] # which patches were included in the result
483 for i, (dataId, maskedImage) in enumerate(maskedImages.items()):
484 # Only merge into the trimmed box, to save memory
485 clippedBox = geom.Box2I(maskedImage.getBBox())
486 clippedBox.clip(bbox)
487 if clippedBox.area == 0:
488 self.log.debug("%s does not overlap template region.", dataId)
489 continue # nothing in this image overlaps the output
490 maskedImage = maskedImage.subset(clippedBox)
491 # Catch both zero-value and NaN variance plane pixels
492 good = (maskedImage.variance.array > 0) & (
493 np.isfinite(maskedImage.variance.array)
494 )
495 weight = maskedImage.variance.array[good] ** (-0.5)
496 bad = np.isnan(maskedImage.image.array) | ~good
497 # Note that modifying the patch MaskedImage in place is fine;
498 # we're throwing it away at the end anyway.
499 maskedImage.image.array[bad] = 0.0
500 maskedImage.variance.array[bad] = 0.0
501 # Reset mask, too, since these pixels don't contribute to sum.
502 maskedImage.mask.array[bad] = 0
503 # Cannot use `merged.maskedImage *= weight` because that operator
504 # multiplies the variance by the weight twice; in this case
505 # `weight` are the exact values we want to scale by.
506 maskedImage.image.array[good] *= weight
507 maskedImage.variance.array[good] *= weight
508 weights[clippedBox].array[good] += weight
509 # Free memory before creating new large arrays
510 del weight
511 merged.maskedImage[clippedBox] += maskedImage
512 included.append(i)
513
514 good = weights.array > 0
515
516 # Cannot use `merged.maskedImage /= weights` because that
517 # operator divides the variance by the weight twice; in this case
518 # `weights` are the exact values we want to scale by.
519 weights = weights.array[good]
520 merged.image.array[good] /= weights
521 merged.variance.array[good] /= weights
522
523 merged.mask.array[~good] |= merged.mask.getPlaneBitMask("NO_DATA")
524
525 return merged, good.sum(), included
526

◆ _selectDataRef()

lsst.ip.diffim.getTemplate._selectDataRef ( coaddRef,
tract,
patch )
protected

Definition at line 805 of file getTemplate.py.

805def _selectDataRef(coaddRef, tract, patch):
806 condition = (coaddRef.dataId["tract"] == tract) & (
807 coaddRef.dataId["patch"] == patch
808 )
809 return condition

◆ checkPatchList()

lsst.ip.diffim.getTemplate.checkPatchList ( self,
patchList )
Check that all of the DcrModel subfilters are present for each
patch.

Parameters
----------
patchList : `dict`
    Dict of the patches containing valid data for each tract.

Raises
------
RuntimeError
    If the number of exposures found for a patch does not match the
    number of subfilters.

Definition at line 744 of file getTemplate.py.

744 def checkPatchList(self, patchList):
745 """Check that all of the DcrModel subfilters are present for each
746 patch.
747
748 Parameters
749 ----------
750 patchList : `dict`
751 Dict of the patches containing valid data for each tract.
752
753 Raises
754 ------
755 RuntimeError
756 If the number of exposures found for a patch does not match the
757 number of subfilters.
758 """
759 for tract in patchList:
760 for patch in set(patchList[tract]):
761 if patchList[tract].count(patch) != self.config.numSubfilters:
762 raise RuntimeError(
763 "Invalid number of DcrModel subfilters found: %d vs %d expected",
764 patchList[tract].count(patch),
765 self.config.numSubfilters,
766 )
767

◆ getDcrModel()

lsst.ip.diffim.getTemplate.getDcrModel ( self,
patchList,
coaddRefs,
visitInfo )
Build DCR-matched coadds from a list of exposure references.

Parameters
----------
patchList : `dict`
    Dict of the patches containing valid data for each tract.
coaddRefs : `list` [`lsst.daf.butler.DeferredDatasetHandle`]
    Data references to `~lsst.afw.image.Exposure` representing
    DcrModels that overlap the detector.
visitInfo : `lsst.afw.image.VisitInfo`
    Metadata for the science image.

Returns
-------
coaddExposures : `list` [`lsst.afw.image.Exposure`]
    Coadd exposures that overlap the detector.

Definition at line 768 of file getTemplate.py.

768 def getDcrModel(self, patchList, coaddRefs, visitInfo):
769 """Build DCR-matched coadds from a list of exposure references.
770
771 Parameters
772 ----------
773 patchList : `dict`
774 Dict of the patches containing valid data for each tract.
775 coaddRefs : `list` [`lsst.daf.butler.DeferredDatasetHandle`]
776 Data references to `~lsst.afw.image.Exposure` representing
777 DcrModels that overlap the detector.
778 visitInfo : `lsst.afw.image.VisitInfo`
779 Metadata for the science image.
780
781 Returns
782 -------
783 coaddExposures : `list` [`lsst.afw.image.Exposure`]
784 Coadd exposures that overlap the detector.
785 """
786 coaddExposures = collections.defaultdict(list)
787 for tract in patchList:
788 for patch in set(patchList[tract]):
789 coaddRefList = [
790 coaddRef
791 for coaddRef in coaddRefs
792 if _selectDataRef(coaddRef, tract, patch)
793 ]
794
795 dcrModel = DcrModel.fromQuantum(
796 coaddRefList,
797 self.config.effectiveWavelength,
798 self.config.bandwidth,
799 self.config.numSubfilters,
800 )
801 coaddExposures[tract].append(dcrModel.buildMatchedExposureHandle(visitInfo=visitInfo))
802 return coaddExposures
803
804

◆ run()

lsst.ip.diffim.getTemplate.run ( self,
* ,
coaddExposureHandles,
bbox,
wcs,
dataIds,
physical_filter )
Warp coadds from multiple tracts and patches to form a template to
subtract from a science image.

Tract and patch overlap regions are combined by a variance-weighted
average, and the variance planes are combined with the same weights,
not added in quadrature; the overlap regions are not statistically
independent, because they're derived from the same original data.
The PSF on the template is created by combining the CoaddPsf on each
template image into a meta-CoaddPsf.

Parameters
----------
coaddExposureHandles : `dict` [`int`,  `list` of \
                  [`lsst.daf.butler.DeferredDatasetHandle` of \
                   `lsst.afw.image.Exposure`]]
    Coadds to be mosaicked, indexed on tract id.
bbox : `lsst.geom.Box2I`
    Template Bounding box of the detector geometry onto which to
    resample the ``coaddExposureHandles``. Modified in-place to include the
    template border.
wcs : `lsst.afw.geom.SkyWcs`
    Template WCS onto which to resample the ``coaddExposureHandles``.
dataIds : `dict` [`int`, `list` [`lsst.daf.butler.DataCoordinate`]]
    Record of the tract and patch of each coaddExposure, indexed on
    tract id.
physical_filter : `str`
    Physical filter of the science image.

Returns
-------
result : `lsst.pipe.base.Struct`
   A struct with attributes:

   ``template``
       A template coadd exposure assembled out of patches
       (`lsst.afw.image.ExposureF`).

Raises
------
NoWorkFound
    If no coadds are found with sufficient un-masked pixels.

Definition at line 255 of file getTemplate.py.

255 def run(self, *, coaddExposureHandles, bbox, wcs, dataIds, physical_filter):
256 """Warp coadds from multiple tracts and patches to form a template to
257 subtract from a science image.
258
259 Tract and patch overlap regions are combined by a variance-weighted
260 average, and the variance planes are combined with the same weights,
261 not added in quadrature; the overlap regions are not statistically
262 independent, because they're derived from the same original data.
263 The PSF on the template is created by combining the CoaddPsf on each
264 template image into a meta-CoaddPsf.
265
266 Parameters
267 ----------
268 coaddExposureHandles : `dict` [`int`, `list` of \
269 [`lsst.daf.butler.DeferredDatasetHandle` of \
270 `lsst.afw.image.Exposure`]]
271 Coadds to be mosaicked, indexed on tract id.
272 bbox : `lsst.geom.Box2I`
273 Template Bounding box of the detector geometry onto which to
274 resample the ``coaddExposureHandles``. Modified in-place to include the
275 template border.
276 wcs : `lsst.afw.geom.SkyWcs`
277 Template WCS onto which to resample the ``coaddExposureHandles``.
278 dataIds : `dict` [`int`, `list` [`lsst.daf.butler.DataCoordinate`]]
279 Record of the tract and patch of each coaddExposure, indexed on
280 tract id.
281 physical_filter : `str`
282 Physical filter of the science image.
283
284 Returns
285 -------
286 result : `lsst.pipe.base.Struct`
287 A struct with attributes:
288
289 ``template``
290 A template coadd exposure assembled out of patches
291 (`lsst.afw.image.ExposureF`).
292
293 Raises
294 ------
295 NoWorkFound
296 If no coadds are found with sufficient un-masked pixels.
297 """
298 band, photoCalib = self._checkInputs(dataIds, coaddExposureHandles)
299
300 bbox.grow(self.config.templateBorderSize)
301
302 warped = {}
303 catalogs = []
304 for tract in coaddExposureHandles:
305 maskedImages, catalog, totalBox = self._makeExposureCatalog(
306 coaddExposureHandles[tract], dataIds[tract]
307 )
308 warpedBox = computeWarpedBBox(catalog[0].wcs, bbox, wcs)
309 warpedBox.grow(5) # to ensure we catch all relevant input pixels
310 # Combine images from individual patches together.
311 unwarped, count, included = self._merge(
312 maskedImages, warpedBox, catalog[0].wcs
313 )
314 # Delete `maskedImages` after combining into one large image to reduce peak memory use
315 del maskedImages
316 if count == 0:
317 self.log.info(
318 "No valid pixels from coadd patches in tract %s; not including in output.",
319 tract,
320 )
321 continue
322 warpedBox.clip(totalBox)
323 potentialInput = self.warper.warpExposure(
324 wcs, unwarped.subset(warpedBox), destBBox=bbox
325 )
326
327 # Delete the single large `unwarped` image after warping to reduce peak memory use
328 del unwarped
329 if np.all(
330 potentialInput.mask.array
331 & potentialInput.mask.getPlaneBitMask("NO_DATA")
332 ):
333 self.log.info(
334 "No overlap from coadd patches in tract %s; not including in output.",
335 tract,
336 )
337 continue
338
339 # Trim the exposure catalog to just the patches that were used.
340 tempCatalog = afwTable.ExposureCatalog(self.schema)
341 tempCatalog.reserve(len(included))
342 for i in included:
343 tempCatalog.append(catalog[i])
344 catalogs.append(tempCatalog)
345 warped[tract] = potentialInput.maskedImage
346
347 if len(warped) == 0:
348 raise pipeBase.NoWorkFound("No patches found to overlap science exposure.")
349 # At this point, all entries will be valid, so we can ignore included.
350 template, count, _ = self._merge(warped, bbox, wcs)
351 if count == 0:
352 raise pipeBase.NoWorkFound("No valid pixels in warped template.")
353
354 # Make a single catalog containing all the inputs that were accepted.
355 catalog = afwTable.ExposureCatalog(self.schema)
356 catalog.reserve(sum([len(c) for c in catalogs]))
357 for c in catalogs:
358 catalog.extend(c)
359
360 template.setPsf(self._makePsf(template, catalog, wcs))
361 template.setFilter(afwImage.FilterLabel(band, physical_filter))
362 template.setPhotoCalib(photoCalib)
363 return pipeBase.Struct(template=template)
364
A group of labels for a filter in an exposure or coadd.
Definition FilterLabel.h:58

Variable Documentation

◆ coaddExposures

lsst.ip.diffim.getTemplate.coaddExposures = collections.defaultdict(list)

Definition at line 228 of file getTemplate.py.

◆ dataId

lsst.ip.diffim.getTemplate.dataId = coaddRef.dataId

Definition at line 232 of file getTemplate.py.

◆ dataIds

lsst.ip.diffim.getTemplate.dataIds = collections.defaultdict(list)

Definition at line 229 of file getTemplate.py.

◆ detectorPolygon

lsst.ip.diffim.getTemplate.detectorPolygon = geom.Box2D(bbox)

Definition at line 226 of file getTemplate.py.

◆ overlappingArea

int lsst.ip.diffim.getTemplate.overlappingArea = 0

Definition at line 227 of file getTemplate.py.

◆ patchBBox

lsst.ip.diffim.getTemplate.patchBBox = skymap[dataId["tract"]][dataId["patch"]].getOuterBBox()

Definition at line 234 of file getTemplate.py.

◆ patchCorners

lsst.ip.diffim.getTemplate.patchCorners = patchWcs.pixelToSky(geom.Box2D(patchBBox).getCorners())

Definition at line 235 of file getTemplate.py.

◆ patchList

lsst.ip.diffim.getTemplate.patchList = dict()

Definition at line 711 of file getTemplate.py.

◆ patchPolygon

lsst.ip.diffim.getTemplate.patchPolygon = afwGeom.Polygon(wcs.skyToPixel(patchCorners))

Definition at line 236 of file getTemplate.py.

◆ patchWcs

lsst.ip.diffim.getTemplate.patchWcs = skymap[dataId["tract"]].getWcs()

Definition at line 233 of file getTemplate.py.

◆ subfilter

lsst.ip.diffim.getTemplate.subfilter = dataId["subfilter"]

Definition at line 714 of file getTemplate.py.