Coverage for python/lsst/pipe/tasks/multiBand.py: 23%

1# This file is part of pipe_tasks. 

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__ = ["DetectCoaddSourcesConfig", "DetectCoaddSourcesTask", 

23 "MeasureMergedCoaddSourcesConfig", "MeasureMergedCoaddSourcesTask", 

24 ] 

25 

26import numpy as np 

27 

28from lsst.geom import Extent2I 

29from lsst.pipe.base import ( 

30 AnnotatedPartialOutputsError, 

31 Struct, 

32 PipelineTask, 

33 PipelineTaskConfig, 

34 PipelineTaskConnections 

35) 

36import lsst.pipe.base.connectionTypes as cT 

37from lsst.pex.config import Field, ConfigurableField, ChoiceField 

38from lsst.meas.algorithms import ( 

39 DynamicDetectionTask, 

40 ExceedsMaxVarianceScaleError, 

41 InsufficientSourcesError, 

42 ReferenceObjectLoader, 

43 ScaleVarianceTask, 

44 SetPrimaryFlagsTask, 

45 TooManyMaskedPixelsError, 

46) 

47from lsst.meas.base import ( 

48 SingleFrameMeasurementTask, 

49 ApplyApCorrTask, 

50 CatalogCalculationTask, 

51 SkyMapIdGeneratorConfig, 

52) 

53from lsst.meas.extensions.scarlet.io import updateCatalogFootprints 

54from lsst.meas.astrom import DirectMatchTask, denormalizeMatches 

55from lsst.pipe.tasks.propagateSourceFlags import PropagateSourceFlagsTask 

56import lsst.afw.image as afwImage 

57import lsst.afw.math as afwMath 

58import lsst.afw.table as afwTable 

59from lsst.daf.base import PropertyList 

60from lsst.skymap import BaseSkyMap 

61 

62# NOTE: these imports are a convenience so multiband users only have to import this file. 

63from .mergeDetections import MergeDetectionsConfig, MergeDetectionsTask # noqa: F401 

64from .mergeMeasurements import MergeMeasurementsConfig, MergeMeasurementsTask # noqa: F401 

65from .multiBandUtils import CullPeaksConfig # noqa: F401 

66from .deblendCoaddSourcesPipeline import DeblendCoaddSourcesMultiConfig # noqa: F401 

67from .deblendCoaddSourcesPipeline import DeblendCoaddSourcesMultiTask # noqa: F401 

68 

69 

70""" 

71New set types: 

72* deepCoadd_det: detections from what used to be processCoadd (tract, patch, filter) 

73* deepCoadd_mergeDet: merged detections (tract, patch) 

74* deepCoadd_meas: measurements of merged detections (tract, patch, filter) 

75* deepCoadd_ref: reference sources (tract, patch) 

76All of these have associated *_schema catalogs that require no data ID and hold no records. 

77 

78In addition, we have a schema-only dataset, which saves the schema for the PeakRecords in 

79the mergeDet, meas, and ref dataset Footprints: 

80* deepCoadd_peak_schema 

81""" 

82 

83 

84############################################################################################################## 

85class DetectCoaddSourcesConnections(PipelineTaskConnections, 

86 dimensions=("tract", "patch", "band", "skymap"), 

87 defaultTemplates={"inputCoaddName": "deep", "outputCoaddName": "deep"}): 

88 detectionSchema = cT.InitOutput( 

89 doc="Schema of the detection catalog", 

90 name="{outputCoaddName}Coadd_det_schema", 

91 storageClass="SourceCatalog", 

92 ) 

93 exposure = cT.Input( 

94 doc="Exposure on which detections are to be performed. ", 

95 name="{inputCoaddName}Coadd", 

96 storageClass="ExposureF", 

97 dimensions=("tract", "patch", "band", "skymap") 

98 ) 

99 exposure_cells = cT.Input( 

100 doc="Exposure on which detections are to be performed. ", 

101 name="{inputCoaddName}CoaddCell", 

102 storageClass="MultipleCellCoadd", 

103 dimensions=("tract", "patch", "band", "skymap"), 

104 ) 

105 skyMap = cT.Input( 

106 doc="Description of the skymap's tracts and patches.", 

107 name=BaseSkyMap.SKYMAP_DATASET_TYPE_NAME, 

108 storageClass="SkyMap", 

109 dimensions=("skymap",), 

110 ) 

111 outputBackgrounds = cT.Output( 

112 doc="Output Backgrounds used in detection", 

113 name="{outputCoaddName}Coadd_calexp_background", 

114 storageClass="Background", 

115 dimensions=("tract", "patch", "band", "skymap") 

116 ) 

117 outputSources = cT.Output( 

118 doc="Detected sources catalog", 

119 name="{outputCoaddName}Coadd_det", 

120 storageClass="SourceCatalog", 

121 dimensions=("tract", "patch", "band", "skymap") 

122 ) 

123 outputExposure = cT.Output( 

124 doc="Exposure post detection", 

125 name="{outputCoaddName}Coadd_calexp", 

126 storageClass="ExposureF", 

127 dimensions=("tract", "patch", "band", "skymap") 

128 ) 

129 

130 def __init__(self, *, config=None): 

131 super().__init__(config=config) 

132 assert isinstance(config, DetectCoaddSourcesConfig) 

133 

134 if config.useCellCoadds: 

135 del self.exposure 

136 else: 

137 del self.exposure_cells 

138 

139 if not self.config.forceExactBinning: 

140 del self.skyMap 

141 if self.config.writeOnlyBackgrounds: 

142 del self.outputExposure 

143 del self.outputSources 

144 del self.detectionSchema 

145 

146 

147class DetectCoaddSourcesConfig(PipelineTaskConfig, pipelineConnections=DetectCoaddSourcesConnections): 

148 """Configuration parameters for the DetectCoaddSourcesTask 

149 """ 

150 

151 doScaleVariance = Field(dtype=bool, default=True, doc="Scale variance plane using empirical noise?") 

152 scaleVariance = ConfigurableField(target=ScaleVarianceTask, doc="Variance rescaling") 

153 detection = ConfigurableField(target=DynamicDetectionTask, doc="Source detection") 

154 coaddName = Field(dtype=str, default="deep", doc="Name of coadd") 

155 useCellCoadds = Field(dtype=bool, default=False, doc="Whether to use cell coadds?") 

156 hasFakes = Field( 

157 dtype=bool, 

158 default=False, 

159 doc="Should be set to True if fake sources have been inserted into the input data.", 

160 ) 

161 idGenerator = SkyMapIdGeneratorConfig.make_field() 

162 forceExactBinning = Field( 

163 dtype=bool, 

164 default=False, 

165 doc=( 

166 "Check that the background bin size evenly divides the patch inner region, and " 

167 "crop the outer region to an integer number of bins." 

168 ) 

169 ) 

170 writeOnlyBackgrounds = Field(dtype=bool, default=False, doc="If true, only save the background models.") 

171 writeEmptyBackgrounds = Field( 

172 dtype=bool, 

173 default=False, 

174 doc=( 

175 "If true, save a placeholder background with NaNs in all bins (but the right geometry) when " 

176 "there are no pixels to compute a background from. This can be useful if a later task combines " 

177 "backgrounds from multiple patches as input." 

178 ) 

179 ) 

180 

181 def setDefaults(self): 

182 super().setDefaults() 

183 self.detection.thresholdType = "pixel_stdev" 

184 self.detection.isotropicGrow = True 

185 # Coadds are made from background-subtracted CCDs, so any background subtraction should be very basic 

186 self.detection.reEstimateBackground = False 

187 self.detection.background.useApprox = False 

188 self.detection.background.binSize = 4096 

189 self.detection.background.undersampleStyle = 'REDUCE_INTERP_ORDER' 

190 self.detection.doTempWideBackground = True # Suppress large footprints that overwhelm the deblender 

191 # Include band in packed data IDs that go into object IDs (None -> "as 

192 # many bands as are defined", rather than the default of zero). 

193 self.idGenerator.packer.n_bands = None 

194 

195 

196class DetectCoaddSourcesTask(PipelineTask): 

197 """Detect sources on a single filter coadd. 

198 

199 Coadding individual visits requires each exposure to be warped. This 

200 introduces covariance in the noise properties across pixels. Before 

201 detection, we correct the coadd variance by scaling the variance plane in 

202 the coadd to match the observed variance. This is an approximate 

203 approach -- strictly, we should propagate the full covariance matrix -- 

204 but it is simple and works well in practice. 

205 

206 After scaling the variance plane, we detect sources and generate footprints 

207 by delegating to the @ref SourceDetectionTask_ "detection" subtask. 

208 

209 DetectCoaddSourcesTask is meant to be run after assembling a coadded image 

210 in a given band. The purpose of the task is to update the background, 

211 detect all sources in a single band and generate a set of parent 

212 footprints. Subsequent tasks in the multi-band processing procedure will 

213 merge sources across bands and, eventually, perform forced photometry. 

214 

215 Parameters 

216 ---------- 

217 schema : `lsst.afw.table.Schema`, optional 

218 Initial schema for the output catalog, modified-in place to include all 

219 fields set by this task. If None, the source minimal schema will be used. 

220 **kwargs 

221 Additional keyword arguments. 

222 """ 

223 

224 _DefaultName = "detectCoaddSources" 

225 ConfigClass = DetectCoaddSourcesConfig 

226 

227 def __init__(self, schema=None, **kwargs): 

228 # N.B. Super is used here to handle the multiple inheritance of PipelineTasks, the init tree 

229 # call structure has been reviewed carefully to be sure super will work as intended. 

230 super().__init__(**kwargs) 

231 if schema is None: 

232 schema = afwTable.SourceTable.makeMinimalSchema() 

233 self.schema = schema 

234 self.makeSubtask("detection", schema=self.schema) 

235 if self.config.doScaleVariance: 

236 self.makeSubtask("scaleVariance") 

237 

238 self.detectionSchema = afwTable.SourceCatalog(self.schema) 

239 

240 def runQuantum(self, butlerQC, inputRefs, outputRefs): 

241 inputs = butlerQC.get(inputRefs) 

242 idGenerator = self.config.idGenerator.apply(butlerQC.quantum.dataId) 

243 

244 if self.config.useCellCoadds: 

245 multiple_cell_coadd = inputs.pop("exposure_cells") 

246 exposure = multiple_cell_coadd.stitch().asExposure() 

247 else: 

248 exposure = inputs.pop("exposure") 

249 

250 skyMap = inputs.pop("skyMap", None) 

251 if skyMap is not None: 

252 patchInfo = skyMap[butlerQC.quantum.dataId["tract"]][butlerQC.quantum.dataId["patch"]] 

253 else: 

254 patchInfo = None 

255 

256 assert not inputs, "runQuantum got more inputs than expected." 

257 try: 

258 outputs = self.run( 

259 exposure=exposure, 

260 idFactory=idGenerator.make_table_id_factory(), 

261 expId=idGenerator.catalog_id, 

262 patchInfo=patchInfo, 

263 ) 

264 except (TooManyMaskedPixelsError, ExceedsMaxVarianceScaleError, InsufficientSourcesError) as e: 

265 if self.config.writeEmptyBackgrounds: 

266 butlerQC.put(self._makeEmptyBackground(exposure, patchInfo), outputRefs.outputBackgrounds) 

267 error = AnnotatedPartialOutputsError.annotate( 

268 e, 

269 self, 

270 exposure, 

271 log=self.log, 

272 ) 

273 raise error from e 

274 

275 butlerQC.put(outputs, outputRefs) 

276 

277 def run(self, exposure, idFactory, expId, patchInfo=None): 

278 """Run detection on an exposure. 

279 

280 First scale the variance plane to match the observed variance 

281 using ``ScaleVarianceTask``. Then invoke the ``SourceDetectionTask_`` "detection" subtask to 

282 detect sources. 

283 

284 Parameters 

285 ---------- 

286 exposure : `lsst.afw.image.Exposure` 

287 Exposure on which to detect (may be background-subtracted and scaled, 

288 depending on configuration). 

289 idFactory : `lsst.afw.table.IdFactory` 

290 IdFactory to set source identifiers. 

291 expId : `int` 

292 Exposure identifier (integer) for RNG seed. 

293 patchInfo : `lsst.skymap.PatchInfo`, optional 

294 Description of the patch geometry. Only needed if 

295 `~DetectCoaddSourceConfig.forceExactBinning` is `True`. 

296 

297 Returns 

298 ------- 

299 result : `lsst.pipe.base.Struct` 

300 Results as a struct with attributes: 

301 

302 ``sources`` 

303 Catalog of detections (`lsst.afw.table.SourceCatalog`). 

304 ``backgrounds`` 

305 List of backgrounds (`list`). 

306 """ 

307 if self.config.forceExactBinning: 

308 exposure = self._cropToExactBinning(exposure, patchInfo) 

309 if self.config.doScaleVariance: 

310 varScale = self.scaleVariance.run(exposure.maskedImage) 

311 exposure.getMetadata().add("VARIANCE_SCALE", varScale) 

312 backgrounds = afwMath.BackgroundList() 

313 table = afwTable.SourceTable.make(self.schema, idFactory) 

314 detections = self.detection.run(table, exposure, expId=expId) 

315 sources = detections.sources 

316 if hasattr(detections, "background") and detections.background: 

317 for bg in detections.background: 

318 backgrounds.append(bg) 

319 if len(backgrounds) == 0: 

320 # Persist a constant background with value of 0.0 to get around 

321 # inability to persist empty BackgroundList. 

322 bgLevel = 0.0 

323 bgStats = afwImage.MaskedImageF(1, 1) 

324 bgStats.set(bgLevel, 0, bgLevel) 

325 bg = afwMath.BackgroundMI(exposure.getBBox(), bgStats) 

326 bgData = (bg, afwMath.Interpolate.LINEAR, afwMath.REDUCE_INTERP_ORDER, 

327 afwMath.ApproximateControl.UNKNOWN, 0, 0, False) 

328 backgrounds.append(bgData) 

329 

330 return Struct(outputSources=sources, outputBackgrounds=backgrounds, outputExposure=exposure) 

331 

332 def _cropToExactBinning(self, exposure, patchInfo): 

333 """Crop a coadd `~lsst.afw.image.Exposure` instance to ensure exact 

334 background binning. 

335 

336 Parameters 

337 ---------- 

338 exposure : `lsst.afw.image.Exposure` 

339 Exposure to crop, assumed to cover the patch outer bounding box. 

340 patchInfo : `lsst.skymap.PatchInfo` 

341 Description of the patch geometry. 

342 

343 Returns 

344 ------- 

345 cropped : `lsst.afw.image.Exposure` 

346 View of ``exposure`` with background bins that evenly divide both 

347 the full cropped image and the patch inner region. The bounding 

348 box is guaranteed to contain the patch inner bounding box and be 

349 contained by the patch outer bounding box. 

350 

351 Raises 

352 ------ 

353 ValueError 

354 Raised if the patch inner region width or height is not a multiple 

355 of the background bin size. 

356 """ 

357 bbox = patchInfo.getInnerBBox() 

358 if bbox.width % self.detection.background.binSizeX: 

359 raise ValueError( 

360 f"Patch inner width {bbox.width} does not evenly " 

361 f"divide bin width {self.detection.background.binSizeX}." 

362 ) 

363 if bbox.height % self.detection.background.binSizeY: 

364 raise ValueError( 

365 f"Patch inner height {bbox.height} does not evenly " 

366 f"divide bin height {self.detection.background.binSizeY}." 

367 ) 

368 outer_bbox = patchInfo.getOuterBBox() 

369 n_bins_grow_x = (bbox.x.begin - outer_bbox.x.begin) // self.detection.background.binSizeX 

370 n_bins_grow_y = (bbox.y.begin - outer_bbox.y.begin) // self.detection.background.binSizeY 

371 bbox.grow( 

372 Extent2I( 

373 n_bins_grow_x*self.detection.background.binSizeX, 

374 n_bins_grow_y*self.detection.background.binSizeY, 

375 ) 

376 ) 

377 assert outer_bbox.contains(bbox) 

378 assert bbox.contains(patchInfo.getInnerBBox()) 

379 assert bbox.width % self.detection.background.binSizeX == 0 

380 assert bbox.height % self.detection.background.binSizeY == 0 

381 return exposure[bbox] 

382 

383 def _makeEmptyBackground(self, exposure, patchInfo=None): 

384 """Construct an empty `lsst.afw.math.BackgroundList` with NaN values. 

385 

386 Parameters 

387 ---------- 

388 exposure : `lsst.afw.image.Exposure` 

389 Exposure that the background should correspond to. 

390 patchInfo : `lsst.skymap.PatchInfo`, optional 

391 Description of the patch geometry. Only needed if 

392 `~DetectCoaddSourceConfig.forceExactBinning` is `True`. 

393 

394 Returns 

395 ------- 

396 background : `lsst.afw.math.BackgroundList` 

397 A background object with a single layer and the same bin geometry 

398 that a background for that exposure would have had if it had enough 

399 usable pixels. This object cannot actually be used for background 

400 subtraction. 

401 """ 

402 # BackgroundList objects are a huge pain to construct without actually 

403 # measuring the background, so that's what we do - on an image with 

404 # 0s everywhere and no masks. And then we replace the background 

405 # "stats image" with NaNs and NO_DATA. 

406 if self.config.forceExactBinning: 

407 exposure = self._cropToExactBinning(exposure, patchInfo).clone() 

408 exposure.image.array[:, :] = 0.0 

409 exposure.mask[:, :] = 0 

410 background = self.detection.background.run(exposure).background 

411 for bg, *_ in background: 

412 stats = bg.getStatsImage() 

413 stats.image.array[:, :] = np.nan 

414 stats.mask.array[:, :] = stats.mask.getPlaneBitMask("NO_DATA") 

415 stats.variance.array[:, :] = 0.0 

416 return background 

417 

418 

419class MeasureMergedCoaddSourcesConnections( 

420 PipelineTaskConnections, 

421 dimensions=("tract", "patch", "band", "skymap"), 

422 defaultTemplates={ 

423 "inputCoaddName": "deep", 

424 "outputCoaddName": "deep", 

425 "deblendedCatalog": "deblendedFlux", 

426 }, 

427 deprecatedTemplates={ 

428 # TODO[DM-47797]: remove this deprecated connection template. 

429 "deblendedCatalog": "Support for old deblender outputs will be removed after v29." 

430 }, 

431): 

432 inputSchema = cT.InitInput( 

433 doc="Input schema for measure merged task produced by a deblender or detection task", 

434 name="{inputCoaddName}Coadd_deblendedFlux_schema", 

435 storageClass="SourceCatalog" 

436 ) 

437 outputSchema = cT.InitOutput( 

438 doc="Output schema after all new fields are added by task", 

439 name="{inputCoaddName}Coadd_meas_schema", 

440 storageClass="SourceCatalog" 

441 ) 

442 # TODO[DM-47797]: remove this deprecated connection. 

443 refCat = cT.PrerequisiteInput( 

444 doc="Reference catalog used to match measured sources against known sources", 

445 name="ref_cat", 

446 storageClass="SimpleCatalog", 

447 dimensions=("skypix",), 

448 deferLoad=True, 

449 multiple=True, 

450 deprecated="Reference matching in measureCoaddSources will be removed after v29.", 

451 ) 

452 exposure = cT.Input( 

453 doc="Input non-cell-based coadd image", 

454 name="{inputCoaddName}Coadd_calexp", 

455 storageClass="ExposureF", 

456 dimensions=("tract", "patch", "band", "skymap") 

457 ) 

458 exposure_cells = cT.Input( 

459 doc="Input cell-based coadd image", 

460 name="{inputCoaddName}CoaddCell", 

461 storageClass="MultipleCellCoadd", 

462 dimensions=("tract", "patch", "band", "skymap"), 

463 ) 

464 background = cT.Input( 

465 doc="Background to subtract from cell-based coadd image", 

466 name="{inputCoaddName}Coadd_calexp_background", 

467 storageClass="Background", 

468 dimensions=("tract", "patch", "band", "skymap") 

469 ) 

470 skyMap = cT.Input( 

471 doc="SkyMap to use in processing", 

472 name=BaseSkyMap.SKYMAP_DATASET_TYPE_NAME, 

473 storageClass="SkyMap", 

474 dimensions=("skymap",), 

475 ) 

476 # TODO[DM-47424]: remove this deprecated connection. 

477 visitCatalogs = cT.Input( 

478 doc="Deprecated and unused.", 

479 name="src", 

480 dimensions=("instrument", "visit", "detector"), 

481 storageClass="SourceCatalog", 

482 multiple=True, 

483 deprecated="Deprecated and unused. Will be removed after v29.", 

484 ) 

485 sourceTableHandles = cT.Input( 

486 doc=("Source tables that are derived from the ``CalibrateTask`` sources. " 

487 "These tables contain astrometry and photometry flags, and optionally " 

488 "PSF flags."), 

489 name="sourceTable_visit", 

490 storageClass="ArrowAstropy", 

491 dimensions=("instrument", "visit"), 

492 multiple=True, 

493 deferLoad=True, 

494 ) 

495 finalizedSourceTableHandles = cT.Input( 

496 doc=("Finalized source tables from ``FinalizeCalibrationTask``. These " 

497 "tables contain PSF flags from the finalized PSF estimation."), 

498 name="finalized_src_table", 

499 storageClass="ArrowAstropy", 

500 dimensions=("instrument", "visit"), 

501 multiple=True, 

502 deferLoad=True, 

503 ) 

504 finalVisitSummaryHandles = cT.Input( 

505 doc="Final visit summary table", 

506 name="finalVisitSummary", 

507 storageClass="ExposureCatalog", 

508 dimensions=("instrument", "visit"), 

509 multiple=True, 

510 deferLoad=True, 

511 ) 

512 # TODO[DM-47797]: remove this deprecated connection. 

513 inputCatalog = cT.Input( 

514 doc=("Name of the input catalog to use." 

515 "If the single band deblender was used this should be 'deblendedFlux." 

516 "If the multi-band deblender was used this should be 'deblendedModel, " 

517 "or deblendedFlux if the multiband deblender was configured to output " 

518 "deblended flux catalogs. If no deblending was performed this should " 

519 "be 'mergeDet'"), 

520 name="{inputCoaddName}Coadd_{deblendedCatalog}", 

521 storageClass="SourceCatalog", 

522 deprecated="Support for old deblender outputs will be removed after v29.", 

523 dimensions=("tract", "patch", "band", "skymap"), 

524 ) 

525 scarletCatalog = cT.Input( 

526 doc="Catalogs produced by multiband deblending", 

527 name="{inputCoaddName}Coadd_deblendedCatalog", 

528 storageClass="SourceCatalog", 

529 dimensions=("tract", "patch", "skymap"), 

530 ) 

531 scarletModels = cT.Input( 

532 doc="Multiband scarlet models produced by the deblender", 

533 name="{inputCoaddName}Coadd_scarletModelData", 

534 storageClass="LsstScarletModelData", 

535 dimensions=("tract", "patch", "skymap"), 

536 ) 

537 outputSources = cT.Output( 

538 doc="Source catalog containing all the measurement information generated in this task", 

539 name="{outputCoaddName}Coadd_meas", 

540 dimensions=("tract", "patch", "band", "skymap"), 

541 storageClass="SourceCatalog", 

542 ) 

543 # TODO[DM-47797]: remove this deprecated connection. 

544 matchResult = cT.Output( 

545 doc="Match catalog produced by configured matcher, optional on doMatchSources", 

546 name="{outputCoaddName}Coadd_measMatch", 

547 dimensions=("tract", "patch", "band", "skymap"), 

548 storageClass="Catalog", 

549 deprecated="Reference matching in measureCoaddSources will be removed after v29.", 

550 ) 

551 # TODO[DM-47797]: remove this deprecated connection. 

552 denormMatches = cT.Output( 

553 doc="Denormalized Match catalog produced by configured matcher, optional on " 

554 "doWriteMatchesDenormalized", 

555 name="{outputCoaddName}Coadd_measMatchFull", 

556 dimensions=("tract", "patch", "band", "skymap"), 

557 storageClass="Catalog", 

558 deprecated="Reference matching in measureCoaddSources will be removed after v29.", 

559 ) 

560 

561 def __init__(self, *, config=None): 

562 super().__init__(config=config) 

563 del self.visitCatalogs 

564 if not config.doPropagateFlags: 

565 del self.sourceTableHandles 

566 del self.finalizedSourceTableHandles 

567 else: 

568 # Check for types of flags required. 

569 if not config.propagateFlags.source_flags: 

570 del self.sourceTableHandles 

571 if not config.propagateFlags.finalized_source_flags: 

572 del self.finalizedSourceTableHandles 

573 # TODO[DM-47797]: only the 'if' block contents here should survive. 

574 if config.inputCatalog == "deblendedCatalog": 

575 del self.inputCatalog 

576 if not config.doAddFootprints: 

577 del self.scarletModels 

578 else: 

579 del self.deblendedCatalog 

580 del self.scarletModels 

581 

582 # TODO[DM-47797]: delete the conditionals below. 

583 if not config.doMatchSources: 

584 del self.refCat 

585 del self.matchResult 

586 

587 if not config.doWriteMatchesDenormalized: 

588 del self.denormMatches 

589 

590 if config.useCellCoadds: 

591 del self.exposure 

592 else: 

593 del self.exposure_cells 

594 del self.background 

595 

596 

597class MeasureMergedCoaddSourcesConfig(PipelineTaskConfig, 

598 pipelineConnections=MeasureMergedCoaddSourcesConnections): 

599 """Configuration parameters for the MeasureMergedCoaddSourcesTask 

600 """ 

601 inputCatalog = ChoiceField( 

602 dtype=str, 

603 default="deblendedCatalog", 

604 allowed={ 

605 "deblendedCatalog": "Output catalog from ScarletDeblendTask", 

606 "deblendedFlux": "Output catalog from SourceDeblendTask", 

607 "mergeDet": "The merged detections before deblending." 

608 }, 

609 doc="The name of the input catalog.", 

610 # TODO[DM-47797]: remove this config option and anything using it. 

611 deprecated="Support for old deblender outputs will be removed after v29.", 

612 ) 

613 doAddFootprints = Field(dtype=bool, 

614 default=True, 

615 doc="Whether or not to add footprints to the input catalog from scarlet models. " 

616 "This should be true whenever using the multi-band deblender, " 

617 "otherwise this should be False.") 

618 doConserveFlux = Field(dtype=bool, default=True, 

619 doc="Whether to use the deblender models as templates to re-distribute the flux " 

620 "from the 'exposure' (True), or to perform measurements on the deblender " 

621 "model footprints.") 

622 doStripFootprints = Field(dtype=bool, default=True, 

623 doc="Whether to strip footprints from the output catalog before " 

624 "saving to disk. " 

625 "This is usually done when using scarlet models to save disk space.") 

626 useCellCoadds = Field(dtype=bool, default=False, doc="Whether to use cell coadds?") 

627 measurement = ConfigurableField(target=SingleFrameMeasurementTask, doc="Source measurement") 

628 setPrimaryFlags = ConfigurableField(target=SetPrimaryFlagsTask, doc="Set flags for primary tract/patch") 

629 doPropagateFlags = Field( 

630 dtype=bool, default=True, 

631 doc="Whether to match sources to CCD catalogs to propagate flags (to e.g. identify PSF stars)" 

632 ) 

633 propagateFlags = ConfigurableField(target=PropagateSourceFlagsTask, doc="Propagate source flags to coadd") 

634 doMatchSources = Field( 

635 dtype=bool, 

636 default=False, 

637 doc="Match sources to reference catalog?", 

638 deprecated="Reference matching in measureCoaddSources will be removed after v29.", 

639 ) 

640 match = ConfigurableField( 

641 target=DirectMatchTask, 

642 doc="Matching to reference catalog", 

643 deprecated="Reference matching in measureCoaddSources will be removed after v29.", 

644 ) 

645 doWriteMatchesDenormalized = Field( 

646 dtype=bool, 

647 default=False, 

648 doc=("Write reference matches in denormalized format? " 

649 "This format uses more disk space, but is more convenient to read."), 

650 deprecated="Reference matching in measureCoaddSources will be removed after v29.", 

651 ) 

652 coaddName = Field(dtype=str, default="deep", doc="Name of coadd") 

653 psfCache = Field(dtype=int, default=100, doc="Size of psfCache") 

654 checkUnitsParseStrict = Field( 

655 doc="Strictness of Astropy unit compatibility check, can be 'raise', 'warn' or 'silent'", 

656 dtype=str, 

657 default="raise", 

658 ) 

659 doApCorr = Field( 

660 dtype=bool, 

661 default=True, 

662 doc="Apply aperture corrections" 

663 ) 

664 applyApCorr = ConfigurableField( 

665 target=ApplyApCorrTask, 

666 doc="Subtask to apply aperture corrections" 

667 ) 

668 doRunCatalogCalculation = Field( 

669 dtype=bool, 

670 default=True, 

671 doc='Run catalogCalculation task' 

672 ) 

673 catalogCalculation = ConfigurableField( 

674 target=CatalogCalculationTask, 

675 doc="Subtask to run catalogCalculation plugins on catalog" 

676 ) 

677 

678 hasFakes = Field( 

679 dtype=bool, 

680 default=False, 

681 doc="Should be set to True if fake sources have been inserted into the input data." 

682 ) 

683 idGenerator = SkyMapIdGeneratorConfig.make_field() 

684 

685 @property 

686 def refObjLoader(self): 

687 return self.match.refObjLoader 

688 

689 def setDefaults(self): 

690 super().setDefaults() 

691 self.measurement.plugins.names |= ['base_InputCount', 

692 'base_Variance', 

693 'base_LocalPhotoCalib', 

694 'base_LocalWcs'] 

695 

696 # TODO: Remove STREAK in DM-44658, streak masking to happen only in 

697 # ip_diffim; if we can propagate the streak mask from diffim, we can 

698 # still set flags with it here. 

699 self.measurement.plugins['base_PixelFlags'].masksFpAnywhere = ['CLIPPED', 'SENSOR_EDGE', 

700 'INEXACT_PSF'] 

701 self.measurement.plugins['base_PixelFlags'].masksFpCenter = ['CLIPPED', 'SENSOR_EDGE', 

702 'INEXACT_PSF'] 

703 

704 def validate(self): 

705 super().validate() 

706 

707 if not self.doMatchSources and self.doWriteMatchesDenormalized: 

708 raise ValueError("Cannot set doWriteMatchesDenormalized if doMatchSources is False.") 

709 

710 

711class MeasureMergedCoaddSourcesTask(PipelineTask): 

712 """Deblend sources from main catalog in each coadd seperately and measure. 

713 

714 Use peaks and footprints from a master catalog to perform deblending and 

715 measurement in each coadd. 

716 

717 Given a master input catalog of sources (peaks and footprints) or deblender 

718 outputs(including a HeavyFootprint in each band), measure each source on 

719 the coadd. Repeating this procedure with the same master catalog across 

720 multiple coadds will generate a consistent set of child sources. 

721 

722 The deblender retains all peaks and deblends any missing peaks (dropouts in 

723 that band) as PSFs. Source properties are measured and the @c is-primary 

724 flag (indicating sources with no children) is set. Visit flags are 

725 propagated to the coadd sources. 

726 

727 Optionally, we can match the coadd sources to an external reference 

728 catalog. 

729 

730 After MeasureMergedCoaddSourcesTask has been run on multiple coadds, we 

731 have a set of per-band catalogs. The next stage in the multi-band 

732 processing procedure will merge these measurements into a suitable catalog 

733 for driving forced photometry. 

734 

735 Parameters 

736 ---------- 

737 schema : ``lsst.afw.table.Schema`, optional 

738 The schema of the merged detection catalog used as input to this one. 

739 peakSchema : ``lsst.afw.table.Schema`, optional 

740 The schema of the PeakRecords in the Footprints in the merged detection catalog. 

741 refObjLoader : `lsst.meas.algorithms.ReferenceObjectLoader`, optional 

742 An instance of ReferenceObjectLoader that supplies an external reference 

743 catalog. May be None if the loader can be constructed from the butler argument or all steps 

744 requiring a reference catalog are disabled. 

745 initInputs : `dict`, optional 

746 Dictionary that can contain a key ``inputSchema`` containing the 

747 input schema. If present will override the value of ``schema``. 

748 **kwargs 

749 Additional keyword arguments. 

750 """ 

751 

752 _DefaultName = "measureCoaddSources" 

753 ConfigClass = MeasureMergedCoaddSourcesConfig 

754 

755 def __init__(self, schema=None, peakSchema=None, refObjLoader=None, initInputs=None, 

756 **kwargs): 

757 super().__init__(**kwargs) 

758 self.deblended = self.config.inputCatalog.startswith("deblended") 

759 self.inputCatalog = "Coadd_" + self.config.inputCatalog 

760 if initInputs is not None: 

761 schema = initInputs['inputSchema'].schema 

762 if schema is None: 

763 raise ValueError("Schema must be defined.") 

764 self.schemaMapper = afwTable.SchemaMapper(schema) 

765 self.schemaMapper.addMinimalSchema(schema) 

766 self.schema = self.schemaMapper.getOutputSchema() 

767 self.algMetadata = PropertyList() 

768 self.makeSubtask("measurement", schema=self.schema, algMetadata=self.algMetadata) 

769 self.makeSubtask("setPrimaryFlags", schema=self.schema) 

770 # TODO[DM-47797]: remove match subtask 

771 if self.config.doMatchSources: 

772 self.makeSubtask("match", refObjLoader=refObjLoader) 

773 if self.config.doPropagateFlags: 

774 self.makeSubtask("propagateFlags", schema=self.schema) 

775 self.schema.checkUnits(parse_strict=self.config.checkUnitsParseStrict) 

776 if self.config.doApCorr: 

777 self.makeSubtask("applyApCorr", schema=self.schema) 

778 if self.config.doRunCatalogCalculation: 

779 self.makeSubtask("catalogCalculation", schema=self.schema) 

780 

781 self.outputSchema = afwTable.SourceCatalog(self.schema) 

782 

783 def runQuantum(self, butlerQC, inputRefs, outputRefs): 

784 inputs = butlerQC.get(inputRefs) 

785 

786 # TODO[DM-47797]: remove this block 

787 if self.config.doMatchSources: 

788 refObjLoader = ReferenceObjectLoader([ref.datasetRef.dataId for ref in inputRefs.refCat], 

789 inputs.pop('refCat'), 

790 name=self.config.connections.refCat, 

791 config=self.config.refObjLoader, 

792 log=self.log) 

793 self.match.setRefObjLoader(refObjLoader) 

794 

795 if self.config.useCellCoadds: 

796 multiple_cell_coadd = inputs.pop("exposure_cells") 

797 stitched_coadd = multiple_cell_coadd.stitch() 

798 exposure = stitched_coadd.asExposure() 

799 background = inputs.pop("background") 

800 exposure.image -= background.getImage() 

801 

802 ccdInputs = stitched_coadd.ccds 

803 apCorrMap = stitched_coadd.ap_corr_map 

804 band = inputRefs.exposure_cells.dataId["band"] 

805 else: 

806 exposure = inputs.pop("exposure") 

807 # Set psfcache 

808 # move this to run after gen2 deprecation 

809 exposure.getPsf().setCacheCapacity(self.config.psfCache) 

810 

811 ccdInputs = exposure.getInfo().getCoaddInputs().ccds 

812 apCorrMap = exposure.getInfo().getApCorrMap() 

813 band = inputRefs.exposure.dataId["band"] 

814 

815 # Get unique integer ID for IdFactory and RNG seeds; only the latter 

816 # should really be used as the IDs all come from the input catalog. 

817 idGenerator = self.config.idGenerator.apply(butlerQC.quantum.dataId) 

818 

819 # Transform inputCatalog 

820 table = afwTable.SourceTable.make(self.schema, idGenerator.make_table_id_factory()) 

821 sources = afwTable.SourceCatalog(table) 

822 # Load the correct input catalog 

823 if "scarletCatalog" in inputs: 

824 inputCatalog = inputs.pop("scarletCatalog") 

825 catalogRef = inputRefs.scarletCatalog 

826 else: 

827 inputCatalog = inputs.pop("inputCatalog") 

828 catalogRef = inputRefs.inputCatalog 

829 sources.extend(inputCatalog, self.schemaMapper) 

830 del inputCatalog 

831 # Add the HeavyFootprints to the deblended sources 

832 if self.config.doAddFootprints: 

833 modelData = inputs.pop('scarletModels') 

834 if self.config.doConserveFlux: 

835 imageForRedistribution = exposure 

836 else: 

837 imageForRedistribution = None 

838 updateCatalogFootprints( 

839 modelData=modelData, 

840 catalog=sources, 

841 band=band, 

842 imageForRedistribution=imageForRedistribution, 

843 removeScarletData=True, 

844 updateFluxColumns=True, 

845 ) 

846 table = sources.getTable() 

847 table.setMetadata(self.algMetadata) # Capture algorithm metadata to write out to the source catalog. 

848 

849 skyMap = inputs.pop('skyMap') 

850 tractNumber = catalogRef.dataId['tract'] 

851 tractInfo = skyMap[tractNumber] 

852 patchInfo = tractInfo.getPatchInfo(catalogRef.dataId['patch']) 

853 skyInfo = Struct( 

854 skyMap=skyMap, 

855 tractInfo=tractInfo, 

856 patchInfo=patchInfo, 

857 wcs=tractInfo.getWcs(), 

858 bbox=patchInfo.getOuterBBox()