31 from lsst.utils.timer
import timeMethod
33 __all__ = [
"BaseSelectImagesTask",
"BaseExposureInfo",
"WcsSelectImagesTask",
"PsfWcsSelectImagesTask",
34 "DatabaseSelectImagesConfig",
"BestSeeingWcsSelectImagesTask",
"BestSeeingSelectVisitsTask",
35 "BestSeeingQuantileSelectVisitsTask"]
39 """Base configuration for subclasses of BaseSelectImagesTask that use a database"""
40 host = pexConfig.Field(
41 doc=
"Database server host name",
44 port = pexConfig.Field(
45 doc=
"Database server port",
48 database = pexConfig.Field(
49 doc=
"Name of database",
52 maxExposures = pexConfig.Field(
53 doc=
"maximum exposures to select; intended for debugging; ignored if None",
60 """Data about a selected exposure
64 """Create exposure information that can be used to generate data references
66 The object has the following fields:
67 - dataId: data ID of exposure (a dict)
68 - coordList: ICRS coordinates of the corners of the exposure (list of lsst.geom.SpherePoint)
69 plus any others items that are desired
71 super(BaseExposureInfo, self).
__init__(dataId=dataId, coordList=coordList)
75 """Base task for selecting images suitable for coaddition
77 ConfigClass = pexConfig.Config
78 _DefaultName =
"selectImages"
81 def run(self, coordList):
82 """Select images suitable for coaddition in a particular region
84 @param[in] coordList: list of coordinates defining region of interest; if None then select all images
85 subclasses may add additional keyword arguments, as required
87 @return a pipeBase Struct containing:
88 - exposureInfoList: a list of exposure information objects (subclasses of BaseExposureInfo),
89 which have at least the following fields:
90 - dataId: data ID dictionary
91 - coordList: ICRS coordinates of the corners of the exposure (list of lsst.geom.SpherePoint)
93 raise NotImplementedError()
95 def _runArgDictFromDataId(self, dataId):
96 """Extract keyword arguments for run (other than coordList) from a data ID
98 @return keyword arguments for run (other than coordList), as a dict
100 raise NotImplementedError()
102 def runDataRef(self, dataRef, coordList, makeDataRefList=True, selectDataList=[]):
103 """Run based on a data reference
105 This delegates to run() and _runArgDictFromDataId() to do the actual
106 selection. In the event that the selectDataList is non-empty, this will
107 be used to further restrict the selection, providing the user with
108 additional control over the selection.
110 @param[in] dataRef: data reference; must contain any extra keys needed by the subclass
111 @param[in] coordList: list of coordinates defining region of interest; if None, search the whole sky
112 @param[in] makeDataRefList: if True, return dataRefList
113 @param[in] selectDataList: List of SelectStruct with dataRefs to consider for selection
114 @return a pipeBase Struct containing:
115 - exposureInfoList: a list of objects derived from ExposureInfo
116 - dataRefList: a list of data references (None if makeDataRefList False)
119 exposureInfoList = self.
runrun(coordList, **runArgDict).exposureInfoList
121 if len(selectDataList) > 0
and len(exposureInfoList) > 0:
123 ccdKeys, ccdValues = _extractKeyValue(exposureInfoList)
124 inKeys, inValues = _extractKeyValue([s.dataRef
for s
in selectDataList], keys=ccdKeys)
125 inValues =
set(inValues)
126 newExposureInfoList = []
127 for info, ccdVal
in zip(exposureInfoList, ccdValues):
128 if ccdVal
in inValues:
129 newExposureInfoList.append(info)
131 self.log.
info(
"De-selecting exposure %s: not in selectDataList", info.dataId)
132 exposureInfoList = newExposureInfoList
135 butler = dataRef.butlerSubset.butler
136 dataRefList = [butler.dataRef(datasetType=
"calexp",
137 dataId=expInfo.dataId,
138 )
for expInfo
in exposureInfoList]
142 return pipeBase.Struct(
143 dataRefList=dataRefList,
144 exposureInfoList=exposureInfoList,
148 def _extractKeyValue(dataList, keys=None):
149 """Extract the keys and values from a list of dataIds
151 The input dataList is a list of objects that have 'dataId' members.
152 This allows it to be used for both a list of data references and a
155 assert len(dataList) > 0
157 keys = sorted(dataList[0].dataId.keys())
160 for data
in dataList:
161 thisKeys =
set(data.dataId.keys())
162 if thisKeys != keySet:
163 raise RuntimeError(
"DataId keys inconsistent: %s vs %s" % (keySet, thisKeys))
164 values.append(tuple(data.dataId[k]
for k
in keys))
169 """A container for data to be passed to the WcsSelectImagesTask"""
172 super(SelectStruct, self).
__init__(dataRef=dataRef, wcs=wcs, bbox=bbox)
176 """Select images using their Wcs
178 We use the "convexHull" method of lsst.sphgeom.ConvexPolygon to define
179 polygons on the celestial sphere, and test the polygon of the
180 patch for overlap with the polygon of the image.
182 We use "convexHull" instead of generating a ConvexPolygon
183 directly because the standard for the inputs to ConvexPolygon
184 are pretty high and we don't want to be responsible for reaching them.
187 def runDataRef(self, dataRef, coordList, makeDataRefList=True, selectDataList=[]):
188 """Select images in the selectDataList that overlap the patch
190 This method is the old entry point for the Gen2 commandline tasks and drivers
191 Will be deprecated in v22.
193 @param dataRef: Data reference for coadd/tempExp (with tract, patch)
194 @param coordList: List of ICRS coordinates (lsst.geom.SpherePoint) specifying boundary of patch
195 @param makeDataRefList: Construct a list of data references?
196 @param selectDataList: List of SelectStruct, to consider for selection
199 exposureInfoList = []
201 patchVertices = [coord.getVector()
for coord
in coordList]
204 for data
in selectDataList:
205 dataRef = data.dataRef
209 imageCorners = self.
getValidImageCornersgetValidImageCorners(imageWcs, imageBox, patchPoly, dataId=
None)
211 dataRefList.append(dataRef)
214 return pipeBase.Struct(
215 dataRefList=dataRefList
if makeDataRefList
else None,
216 exposureInfoList=exposureInfoList,
219 def run(self, wcsList, bboxList, coordList, dataIds=None, **kwargs):
220 """Return indices of provided lists that meet the selection criteria
224 wcsList : `list` of `lsst.afw.geom.SkyWcs`
225 specifying the WCS's of the input ccds to be selected
226 bboxList : `list` of `lsst.geom.Box2I`
227 specifying the bounding boxes of the input ccds to be selected
228 coordList : `list` of `lsst.geom.SpherePoint`
229 ICRS coordinates specifying boundary of the patch.
233 result: `list` of `int`
234 of indices of selected ccds
237 dataIds = [
None] * len(wcsList)
238 patchVertices = [coord.getVector()
for coord
in coordList]
241 for i, (imageWcs, imageBox, dataId)
in enumerate(zip(wcsList, bboxList, dataIds)):
242 imageCorners = self.
getValidImageCornersgetValidImageCorners(imageWcs, imageBox, patchPoly, dataId)
248 "Return corners or None if bad"
250 imageCorners = [imageWcs.pixelToSky(pix)
for pix
in geom.Box2D(imageBox).getCorners()]
253 self.log.
debug(
"WCS error in testing calexp %s (%s): deselecting", dataId, e)
257 if imagePoly
is None:
258 self.log.
debug(
"Unable to create polygon from image %s: deselecting", dataId)
261 if patchPoly.intersects(imagePoly):
263 self.log.
info(
"Selecting calexp %s", dataId)
268 "Return median absolute deviation scaled to normally distributed data"
269 return 1.4826*np.median(np.abs(array - np.median(array)))
273 dimensions=(
"tract",
"patch",
"skymap",
"instrument",
"visit"),
274 defaultTemplates={
"coaddName":
"deep"}):
278 class PsfWcsSelectImagesConfig(pipeBase.PipelineTaskConfig,
279 pipelineConnections=PsfWcsSelectImagesConnections):
280 maxEllipResidual = pexConfig.Field(
281 doc=
"Maximum median ellipticity residual",
286 maxSizeScatter = pexConfig.Field(
287 doc=
"Maximum scatter in the size residuals",
291 maxScaledSizeScatter = pexConfig.Field(
292 doc=
"Maximum scatter in the size residuals, scaled by the median size",
297 starSelection = pexConfig.Field(
298 doc=
"select star with this field",
300 default=
'calib_psf_used',
301 deprecated=(
'This field has been moved to ComputeExposureSummaryStatsTask and '
302 'will be removed after v24.')
304 starShape = pexConfig.Field(
305 doc=
"name of star shape",
307 default=
'base_SdssShape',
308 deprecated=(
'This field has been moved to ComputeExposureSummaryStatsTask and '
309 'will be removed after v24.')
311 psfShape = pexConfig.Field(
312 doc=
"name of psf shape",
314 default=
'base_SdssShape_psf',
315 deprecated=(
'This field has been moved to ComputeExposureSummaryStatsTask and '
316 'will be removed after v24.')
321 """Select images using their Wcs and cuts on the PSF properties
323 The PSF quality criteria are based on the size and ellipticity residuals from the
324 adaptive second moments of the star and the PSF.
327 - the median of the ellipticty residuals
328 - the robust scatter of the size residuals (using the median absolute deviation)
329 - the robust scatter of the size residuals scaled by the square of
333 ConfigClass = PsfWcsSelectImagesConfig
334 _DefaultName =
"PsfWcsSelectImages"
336 def runDataRef(self, dataRef, coordList, makeDataRefList=True, selectDataList=[]):
337 """Select images in the selectDataList that overlap the patch and satisfy PSF quality critera.
339 This method is the old entry point for the Gen2 commandline tasks and drivers
340 Will be deprecated in v22.
342 @param dataRef: Data reference for coadd/tempExp (with tract, patch)
343 @param coordList: List of ICRS coordinates (lsst.geom.SpherePoint) specifying boundary of patch
344 @param makeDataRefList: Construct a list of data references?
345 @param selectDataList: List of SelectStruct, to consider for selection
347 result = super(PsfWcsSelectImagesTask, self).runDataRef(dataRef, coordList, makeDataRefList,
351 exposureInfoList = []
352 for dataRef, exposureInfo
in zip(result.dataRefList, result.exposureInfoList):
353 butler = dataRef.butlerSubset.butler
354 srcCatalog = butler.get(
'src', dataRef.dataId)
355 valid = self.isValidGen2(srcCatalog, dataRef.dataId)
359 dataRefList.append(dataRef)
360 exposureInfoList.append(exposureInfo)
362 return pipeBase.Struct(
363 dataRefList=dataRefList,
364 exposureInfoList=exposureInfoList,
367 def run(self, wcsList, bboxList, coordList, visitSummary, dataIds=None, **kwargs):
368 """Return indices of provided lists that meet the selection criteria
372 wcsList : `list` of `lsst.afw.geom.SkyWcs`
373 specifying the WCS's of the input ccds to be selected
374 bboxList : `list` of `lsst.geom.Box2I`
375 specifying the bounding boxes of the input ccds to be selected
376 coordList : `list` of `lsst.geom.SpherePoint`
377 ICRS coordinates specifying boundary of the patch.
378 visitSummary : `list` of `lsst.afw.table.ExposureCatalog`
379 containing the PSF shape information for the input ccds to be selected
383 goodPsf: `list` of `int`
384 of indices of selected ccds
386 goodWcs = super(PsfWcsSelectImagesTask, self).
run(wcsList=wcsList, bboxList=bboxList,
387 coordList=coordList, dataIds=dataIds)
391 for i, dataId
in enumerate(dataIds):
394 if self.isValid(visitSummary, dataId[
'detector']):
399 def isValid(self, visitSummary, detectorId):
400 """Should this ccd be selected based on its PSF shape information.
404 visitSummary : `lsst.afw.table.ExposureCatalog`
413 row = visitSummary.find(detectorId)
416 self.log.
warning(
"Removing visit %d detector %d because summary stats not available.",
417 row[
"visit"], detectorId)
420 medianE = np.sqrt(row[
"psfStarDeltaE1Median"]**2. + row[
"psfStarDeltaE2Median"]**2.)
421 scatterSize = row[
"psfStarDeltaSizeScatter"]
422 scaledScatterSize = row[
"psfStarScaledDeltaSizeScatter"]
425 if self.config.maxEllipResidual
and medianE > self.config.maxEllipResidual:
426 self.log.
info(
"Removing visit %d detector %d because median e residual too large: %f vs %f",
427 row[
"visit"], detectorId, medianE, self.config.maxEllipResidual)
429 elif self.config.maxSizeScatter
and scatterSize > self.config.maxSizeScatter:
430 self.log.
info(
"Removing visit %d detector %d because size scatter too large: %f vs %f",
431 row[
"visit"], detectorId, scatterSize, self.config.maxSizeScatter)
433 elif self.config.maxScaledSizeScatter
and scaledScatterSize > self.config.maxScaledSizeScatter:
434 self.log.
info(
"Removing visit %d detector %d because scaled size scatter too large: %f vs %f",
435 row[
"visit"], detectorId, scaledScatterSize, self.config.maxScaledSizeScatter)
440 def isValidGen2(self, srcCatalog, dataId=None):
441 """Should this ccd be selected based on its PSF shape information.
443 This routine is only used in Gen2 processing, and can be
444 removed when Gen2 is retired.
448 srcCatalog : `lsst.afw.table.SourceCatalog`
449 dataId : `dict` of dataId keys, optional.
450 Used only for logging. Defaults to None.
457 mask = srcCatalog[self.config.starSelection]
459 starXX = srcCatalog[self.config.starShape+
'_xx'][mask]
460 starYY = srcCatalog[self.config.starShape+
'_yy'][mask]
461 starXY = srcCatalog[self.config.starShape+
'_xy'][mask]
462 psfXX = srcCatalog[self.config.psfShape+
'_xx'][mask]
463 psfYY = srcCatalog[self.config.psfShape+
'_yy'][mask]
464 psfXY = srcCatalog[self.config.psfShape+
'_xy'][mask]
466 starSize = np.power(starXX*starYY - starXY**2, 0.25)
467 starE1 = (starXX - starYY)/(starXX + starYY)
468 starE2 = 2*starXY/(starXX + starYY)
469 medianSize = np.median(starSize)
471 psfSize = np.power(psfXX*psfYY - psfXY**2, 0.25)
472 psfE1 = (psfXX - psfYY)/(psfXX + psfYY)
473 psfE2 = 2*psfXY/(psfXX + psfYY)
475 medianE1 = np.abs(np.median(starE1 - psfE1))
476 medianE2 = np.abs(np.median(starE2 - psfE2))
477 medianE = np.sqrt(medianE1**2 + medianE2**2)
479 scatterSize =
sigmaMad(starSize - psfSize)
480 scaledScatterSize = scatterSize/medianSize**2
483 if self.config.maxEllipResidual
and medianE > self.config.maxEllipResidual:
484 self.log.
info(
"Removing visit %s because median e residual too large: %f vs %f",
485 dataId, medianE, self.config.maxEllipResidual)
487 elif self.config.maxSizeScatter
and scatterSize > self.config.maxSizeScatter:
488 self.log.
info(
"Removing visit %s because size scatter is too large: %f vs %f",
489 dataId, scatterSize, self.config.maxSizeScatter)
491 elif self.config.maxScaledSizeScatter
and scaledScatterSize > self.config.maxScaledSizeScatter:
492 self.log.
info(
"Removing visit %s because scaled size scatter is too large: %f vs %f",
493 dataId, scaledScatterSize, self.config.maxScaledSizeScatter)
499 class BestSeeingWcsSelectImageConfig(WcsSelectImagesTask.ConfigClass):
500 """Base configuration for BestSeeingSelectImagesTask.
502 nImagesMax = pexConfig.RangeField(
504 doc=
"Maximum number of images to select",
507 maxPsfFwhm = pexConfig.Field(
509 doc=
"Maximum PSF FWHM (in arcseconds) to select",
512 minPsfFwhm = pexConfig.Field(
514 doc=
"Minimum PSF FWHM (in arcseconds) to select",
520 """Select up to a maximum number of the best-seeing images using their Wcs.
522 ConfigClass = BestSeeingWcsSelectImageConfig
524 def runDataRef(self, dataRef, coordList, makeDataRefList=True,
525 selectDataList=None):
526 """Select the best-seeing images in the selectDataList that overlap the patch.
528 This method is the old entry point for the Gen2 commandline tasks and drivers
529 Will be deprecated in v22.
533 dataRef : `lsst.daf.persistence.ButlerDataRef`
534 Data reference for coadd/tempExp (with tract, patch)
535 coordList : `list` of `lsst.geom.SpherePoint`
536 List of ICRS sky coordinates specifying boundary of patch
537 makeDataRefList : `boolean`, optional
538 Construct a list of data references?
539 selectDataList : `list` of `SelectStruct`
540 List of SelectStruct, to consider for selection
544 result : `lsst.pipe.base.Struct`
545 Result struct with components:
546 - ``exposureList``: the selected exposures
547 (`list` of `lsst.pipe.tasks.selectImages.BaseExposureInfo`).
548 - ``dataRefList``: the optional data references corresponding to
549 each element of ``exposureList``
550 (`list` of `lsst.daf.persistence.ButlerDataRef`, or `None`).
554 exposureInfoList = []
556 if selectDataList
is None:
559 result = super().runDataRef(dataRef, coordList, makeDataRefList=
True, selectDataList=selectDataList)
561 for dataRef, exposureInfo
in zip(result.dataRefList, result.exposureInfoList):
562 cal = dataRef.get(
"calexp", immediate=
True)
565 pixToArcseconds = cal.getWcs().getPixelScale().asArcseconds()
566 psfSize = cal.getPsf().computeShape().getDeterminantRadius()*pixToArcseconds
567 sizeFwhm = psfSize * np.sqrt(8.*np.log(2.))
568 if self.config.maxPsfFwhm
and sizeFwhm > self.config.maxPsfFwhm:
570 if self.config.minPsfFwhm
and sizeFwhm < self.config.minPsfFwhm:
572 psfSizes.append(sizeFwhm)
573 dataRefList.append(dataRef)
574 exposureInfoList.append(exposureInfo)
576 if len(psfSizes) > self.config.nImagesMax:
577 sortedIndices = np.argsort(psfSizes)[:self.config.nImagesMax]
578 filteredDataRefList = [dataRefList[i]
for i
in sortedIndices]
579 filteredExposureInfoList = [exposureInfoList[i]
for i
in sortedIndices]
580 self.log.
info(
"%d images selected with FWHM range of %f--%f arcseconds",
581 len(sortedIndices), psfSizes[sortedIndices[0]], psfSizes[sortedIndices[-1]])
584 if len(psfSizes) == 0:
585 self.log.
warning(
"0 images selected.")
587 self.log.
debug(
"%d images selected with FWHM range of %d--%d arcseconds",
588 len(psfSizes), psfSizes[0], psfSizes[-1])
589 filteredDataRefList = dataRefList
590 filteredExposureInfoList = exposureInfoList
592 return pipeBase.Struct(
593 dataRefList=filteredDataRefList
if makeDataRefList
else None,
594 exposureInfoList=filteredExposureInfoList,
598 class BestSeeingSelectVisitsConnections(pipeBase.PipelineTaskConnections,
599 dimensions=(
"tract",
"patch",
"skymap",
"band",
"instrument"),
600 defaultTemplates={
"coaddName":
"goodSeeing"}):
601 skyMap = pipeBase.connectionTypes.Input(
602 doc=
"Input definition of geometry/bbox and projection/wcs for coadded exposures",
603 name=BaseSkyMap.SKYMAP_DATASET_TYPE_NAME,
604 storageClass=
"SkyMap",
605 dimensions=(
"skymap",),
607 visitSummaries = pipeBase.connectionTypes.Input(
608 doc=
"Per-visit consolidated exposure metadata from ConsolidateVisitSummaryTask",
610 storageClass=
"ExposureCatalog",
611 dimensions=(
"instrument",
"visit",),
615 goodVisits = pipeBase.connectionTypes.Output(
616 doc=
"Selected visits to be coadded.",
617 name=
"{coaddName}Visits",
618 storageClass=
"StructuredDataDict",
619 dimensions=(
"instrument",
"tract",
"patch",
"skymap",
"band"),
623 class BestSeeingSelectVisitsConfig(pipeBase.PipelineTaskConfig,
624 pipelineConnections=BestSeeingSelectVisitsConnections):
625 nVisitsMax = pexConfig.RangeField(
627 doc=
"Maximum number of visits to select",
631 maxPsfFwhm = pexConfig.Field(
633 doc=
"Maximum PSF FWHM (in arcseconds) to select",
637 minPsfFwhm = pexConfig.Field(
639 doc=
"Minimum PSF FWHM (in arcseconds) to select",
643 doConfirmOverlap = pexConfig.Field(
645 doc=
"Do remove visits that do not actually overlap the patch?",
648 minMJD = pexConfig.Field(
650 doc=
"Minimum visit MJD to select",
654 maxMJD = pexConfig.Field(
656 doc=
"Maximum visit MJD to select",
662 class BestSeeingSelectVisitsTask(pipeBase.PipelineTask):
663 """Select up to a maximum number of the best-seeing visits
665 Don't exceed the FWHM range specified by configs min(max)PsfFwhm.
666 This Task is a port of the Gen2 image-selector used in the AP pipeline:
667 BestSeeingSelectImagesTask. This Task selects full visits based on the
668 average PSF of the entire visit.
670 ConfigClass = BestSeeingSelectVisitsConfig
671 _DefaultName =
'bestSeeingSelectVisits'
673 def runQuantum(self, butlerQC, inputRefs, outputRefs):
674 inputs = butlerQC.get(inputRefs)
675 quantumDataId = butlerQC.quantum.dataId
676 outputs = self.run(**inputs, dataId=quantumDataId)
677 butlerQC.put(outputs, outputRefs)
679 def run(self, visitSummaries, skyMap, dataId):
684 visitSummary : `list`
685 List of `lsst.pipe.base.connections.DeferredDatasetRef` of
686 visitSummary tables of type `lsst.afw.table.ExposureCatalog`
687 skyMap : `lsst.skyMap.SkyMap`
688 SkyMap for checking visits overlap patch
689 dataId : `dict` of dataId keys
690 For retrieving patch info for checking visits overlap patch
694 result : `lsst.pipe.base.Struct`
695 Result struct with components:
697 - `goodVisits`: `dict` with selected visit ids as keys,
698 so that it can be be saved as a StructuredDataDict.
699 StructuredDataList's are currently limited.
702 if self.config.doConfirmOverlap:
703 patchPolygon = self.makePatchPolygon(skyMap, dataId)
705 inputVisits = [visitSummary.ref.dataId[
'visit']
for visitSummary
in visitSummaries]
708 for visit, visitSummary
in zip(inputVisits, visitSummaries):
710 visitSummary = visitSummary.get()
713 mjd = visitSummary[0].getVisitInfo().getDate().get(system=DateTime.MJD)
715 pixToArcseconds = [vs.getWcs().getPixelScale(vs.getBBox().getCenter()).asArcseconds()
716 for vs
in visitSummary]
718 psfSigmas = np.array([vs[
'psfSigma']
for vs
in visitSummary])
719 fwhm = np.nanmean(psfSigmas * pixToArcseconds) * np.sqrt(8.*np.log(2.))
721 if self.config.maxPsfFwhm
and fwhm > self.config.maxPsfFwhm:
723 if self.config.minPsfFwhm
and fwhm < self.config.minPsfFwhm:
725 if self.config.minMJD
and mjd < self.config.minMJD:
726 self.log.
debug(
'MJD %f earlier than %.2f; rejecting', mjd, self.config.minMJD)
728 if self.config.maxMJD
and mjd > self.config.maxMJD:
729 self.log.
debug(
'MJD %f later than %.2f; rejecting', mjd, self.config.maxMJD)
731 if self.config.doConfirmOverlap
and not self.doesIntersectPolygon(visitSummary, patchPolygon):
734 fwhmSizes.append(fwhm)
737 sortedVisits = [ind
for (_, ind)
in sorted(zip(fwhmSizes, visits))]
738 output = sortedVisits[:self.config.nVisitsMax]
739 self.log.
info(
"%d images selected with FWHM range of %d--%d arcseconds",
740 len(output), fwhmSizes[visits.index(output[0])], fwhmSizes[visits.index(output[-1])])
743 goodVisits = {key:
True for key
in output}
744 return pipeBase.Struct(goodVisits=goodVisits)
746 def makePatchPolygon(self, skyMap, dataId):
747 """Return True if sky polygon overlaps visit
751 skyMap : `lsst.afw.table.ExposureCatalog`
752 Exposure catalog with per-detector geometry
753 dataId : `dict` of dataId keys
754 For retrieving patch info
758 result :` lsst.sphgeom.ConvexPolygon.convexHull`
759 Polygon of patch's outer bbox
761 wcs = skyMap[dataId[
'tract']].getWcs()
762 bbox = skyMap[dataId[
'tract']][dataId[
'patch']].getOuterBBox()
767 def doesIntersectPolygon(self, visitSummary, polygon):
768 """Return True if sky polygon overlaps visit
772 visitSummary : `lsst.afw.table.ExposureCatalog`
773 Exposure catalog with per-detector geometry
774 polygon :` lsst.sphgeom.ConvexPolygon.convexHull`
775 Polygon to check overlap
779 doesIntersect: `bool`
780 Does the visit overlap the polygon
782 doesIntersect =
False
783 for detectorSummary
in visitSummary:
785 zip(detectorSummary[
'raCorners'], detectorSummary[
'decCorners'])]
787 if detectorPolygon.intersects(polygon):
793 class BestSeeingQuantileSelectVisitsConfig(pipeBase.PipelineTaskConfig,
794 pipelineConnections=BestSeeingSelectVisitsConnections):
795 qMin = pexConfig.RangeField(
796 doc=
"Lower bound of quantile range to select. Sorts visits by seeing from narrow to wide, "
797 "and select those in the interquantile range (qMin, qMax). Set qMin to 0 for Best Seeing. "
798 "This config should be changed from zero only for exploratory diffIm testing.",
804 qMax = pexConfig.RangeField(
805 doc=
"Upper bound of quantile range to select. Sorts visits by seeing from narrow to wide, "
806 "and select those in the interquantile range (qMin, qMax). Set qMax to 1 for Worst Seeing.",
812 nVisitsMin = pexConfig.Field(
813 doc=
"At least this number of visits selected and supercedes quantile. For example, if 10 visits "
814 "cover this patch, qMin=0.33, and nVisitsMin=5, the best 5 visits will be selected.",
818 doConfirmOverlap = pexConfig.Field(
820 doc=
"Do remove visits that do not actually overlap the patch?",
823 minMJD = pexConfig.Field(
825 doc=
"Minimum visit MJD to select",
829 maxMJD = pexConfig.Field(
831 doc=
"Maximum visit MJD to select",
837 class BestSeeingQuantileSelectVisitsTask(BestSeeingSelectVisitsTask):
838 """Select a quantile of the best-seeing visits
840 Selects the best (for example, third) full visits based on the average
841 PSF width in the entire visit. It can also be used for difference imaging
842 experiments that require templates with the worst seeing visits.
843 For example, selecting the worst third can be acheived by
844 changing the config parameters qMin to 0.66 and qMax to 1.
846 ConfigClass = BestSeeingQuantileSelectVisitsConfig
847 _DefaultName =
'bestSeeingQuantileSelectVisits'
849 @utils.inheritDoc(BestSeeingSelectVisitsTask)
850 def run(self, visitSummaries, skyMap, dataId):
851 if self.config.doConfirmOverlap:
852 patchPolygon = self.makePatchPolygon(skyMap, dataId)
853 visits = np.array([visitSummary.ref.dataId[
'visit']
for visitSummary
in visitSummaries])
854 radius = np.empty(len(visits))
855 intersects = np.full(len(visits),
True)
856 for i, visitSummary
in enumerate(visitSummaries):
858 visitSummary = visitSummary.get()
860 psfSigma = np.nanmedian([vs[
'psfSigma']
for vs
in visitSummary])
862 if self.config.doConfirmOverlap:
863 intersects[i] = self.doesIntersectPolygon(visitSummary, patchPolygon)
864 if self.config.minMJD
or self.config.maxMJD:
866 mjd = visitSummary[0].getVisitInfo().getDate().get(system=DateTime.MJD)
867 aboveMin = mjd > self.config.minMJD
if self.config.minMJD
else True
868 belowMax = mjd < self.config.maxMJD
if self.config.maxMJD
else True
869 intersects[i] = intersects[i]
and aboveMin
and belowMax
871 sortedVisits = [v
for rad, v
in sorted(zip(radius[intersects], visits[intersects]))]
872 lowerBound =
min(int(np.round(self.config.qMin*len(visits[intersects]))),
873 max(0, len(visits[intersects]) - self.config.nVisitsMin))
874 upperBound =
max(int(np.round(self.config.qMax*len(visits[intersects]))), self.config.nVisitsMin)
877 goodVisits = {int(visit):
True for visit
in sortedVisits[lowerBound:upperBound]}
878 return pipeBase.Struct(goodVisits=goodVisits)
A floating-point coordinate rectangle geometry.
Point in an unspecified spherical coordinate system.
Reports arguments outside the domain of an operation.
Reports errors that are due to events beyond the control of the program.
def __init__(self, dataId, coordList)
def _runArgDictFromDataId(self, dataId)
def runDataRef(self, dataRef, coordList, makeDataRefList=True, selectDataList=[])
def __init__(self, dataRef, wcs, bbox)
def getValidImageCorners(self, imageWcs, imageBox, patchPoly, dataId=None)
def runDataRef(self, dataRef, coordList, makeDataRefList=True, selectDataList=[])
def run(self, wcsList, bboxList, coordList, dataIds=None, **kwargs)
static ConvexPolygon convexHull(std::vector< UnitVector3d > const &points)
convexHull returns the convex hull of the given set of points if it exists and throws an exception ot...
daf::base::PropertyList * list
daf::base::PropertySet * set
def run(self, coaddExposures, bbox, wcs)