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.table as afwTable 

57import lsst.afw.math as afwMath 

58from lsst.daf.base import PropertyList 

59from lsst.skymap import BaseSkyMap 

60 

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

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

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

64from .multiBandUtils import CullPeaksConfig # noqa: F401 

65from .deblendCoaddSourcesPipeline import DeblendCoaddSourcesMultiConfig # noqa: F401 

66from .deblendCoaddSourcesPipeline import DeblendCoaddSourcesMultiTask # noqa: F401 

67 

68 

69""" 

70New set types: 

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

72* deepCoadd_mergeDet: merged detections (tract, patch) 

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

74* deepCoadd_ref: reference sources (tract, patch) 

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

76 

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

78the mergeDet, meas, and ref dataset Footprints: 

79* deepCoadd_peak_schema 

80""" 

81 

82 

83############################################################################################################## 

84class DetectCoaddSourcesConnections(PipelineTaskConnections, 

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

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

87 detectionSchema = cT.InitOutput( 

88 doc="Schema of the detection catalog", 

89 name="{outputCoaddName}Coadd_det_schema", 

90 storageClass="SourceCatalog", 

91 ) 

92 exposure = cT.Input( 

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

94 name="{inputCoaddName}Coadd", 

95 storageClass="ExposureF", 

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

97 ) 

98 exposure_cells = cT.Input( 

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

100 name="{inputCoaddName}CoaddCell", 

101 storageClass="MultipleCellCoadd", 

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

103 ) 

104 skyMap = cT.Input( 

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

106 name=BaseSkyMap.SKYMAP_DATASET_TYPE_NAME, 

107 storageClass="SkyMap", 

108 dimensions=("skymap",), 

109 ) 

110 outputBackgrounds = cT.Output( 

111 doc="Output Backgrounds used in detection", 

112 name="{outputCoaddName}Coadd_calexp_background", 

113 storageClass="Background", 

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

115 ) 

116 outputSources = cT.Output( 

117 doc="Detected sources catalog", 

118 name="{outputCoaddName}Coadd_det", 

119 storageClass="SourceCatalog", 

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

121 ) 

122 outputExposure = cT.Output( 

123 doc="Exposure post detection", 

124 name="{outputCoaddName}Coadd_calexp", 

125 storageClass="ExposureF", 

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

127 ) 

128 

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

130 super().__init__(config=config) 

131 assert isinstance(config, DetectCoaddSourcesConfig) 

132 

133 if config.useCellCoadds: 

134 del self.exposure 

135 else: 

136 del self.exposure_cells 

137 

138 if not self.config.forceExactBinning: 

139 del self.skyMap 

140 if self.config.writeOnlyBackgrounds: 

141 del self.outputExposure 

142 del self.outputSources 

143 del self.detectionSchema 

144 

145 

146class DetectCoaddSourcesConfig(PipelineTaskConfig, pipelineConnections=DetectCoaddSourcesConnections): 

147 """Configuration parameters for the DetectCoaddSourcesTask 

148 """ 

149 

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

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

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

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

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

155 hasFakes = Field( 

156 dtype=bool, 

157 default=False, 

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

159 ) 

160 idGenerator = SkyMapIdGeneratorConfig.make_field() 

161 forceExactBinning = Field( 

162 dtype=bool, 

163 default=False, 

164 doc=( 

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

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

167 ) 

168 ) 

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

170 writeEmptyBackgrounds = Field( 

171 dtype=bool, 

172 default=False, 

173 doc=( 

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

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

176 "backgrounds from multiple patches as input." 

177 ) 

178 ) 

179 

180 def setDefaults(self): 

181 super().setDefaults() 

182 self.detection.thresholdType = "pixel_stdev" 

183 self.detection.isotropicGrow = True 

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

185 self.detection.reEstimateBackground = False 

186 self.detection.background.useApprox = False 

187 self.detection.background.binSize = 4096 

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

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

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

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

192 self.idGenerator.packer.n_bands = None 

193 

194 

195class DetectCoaddSourcesTask(PipelineTask): 

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

197 

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

199 introduces covariance in the noise properties across pixels. Before 

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

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

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

203 but it is simple and works well in practice. 

204 

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

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

207 

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

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

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

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

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

213 

214 Parameters 

215 ---------- 

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

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

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

219 **kwargs 

220 Additional keyword arguments. 

221 """ 

222 

223 _DefaultName = "detectCoaddSources" 

224 ConfigClass = DetectCoaddSourcesConfig 

225 

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

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

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

229 super().__init__(**kwargs) 

230 if schema is None: 

231 schema = afwTable.SourceTable.makeMinimalSchema() 

232 self.schema = schema 

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

234 if self.config.doScaleVariance: 

235 self.makeSubtask("scaleVariance") 

236 

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

238 

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

240 inputs = butlerQC.get(inputRefs) 

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

242 

243 if self.config.useCellCoadds: 

244 multiple_cell_coadd = inputs.pop("exposure_cells") 

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

246 else: 

247 exposure = inputs.pop("exposure") 

248 

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

250 if skyMap is not None: 

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

252 else: 

253 patchInfo = None 

254 

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

256 try: 

257 outputs = self.run( 

258 exposure=exposure, 

259 idFactory=idGenerator.make_table_id_factory(), 

260 expId=idGenerator.catalog_id, 

261 patchInfo=patchInfo, 

262 ) 

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

264 if self.config.writeEmptyBackgrounds: 

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

266 error = AnnotatedPartialOutputsError.annotate( 

267 e, 

268 self, 

269 exposure, 

270 log=self.log, 

271 ) 

272 raise error from e 

273 

274 butlerQC.put(outputs, outputRefs) 

275 

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

277 """Run detection on an exposure. 

278 

279 First scale the variance plane to match the observed variance 

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

281 detect sources. 

282 

283 Parameters 

284 ---------- 

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

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

287 depending on configuration). 

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

289 IdFactory to set source identifiers. 

290 expId : `int` 

291 Exposure identifier (integer) for RNG seed. 

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

293 Description of the patch geometry. Only needed if 

294 `~DetectCoaddSourceConfig.forceExactBinning` is `True`. 

295 

296 Returns 

297 ------- 

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

299 Results as a struct with attributes: 

300 

301 ``sources`` 

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

303 ``backgrounds`` 

304 List of backgrounds (`list`). 

305 """ 

306 if self.config.forceExactBinning: 

307 exposure = self._cropToExactBinning(exposure, patchInfo) 

308 if self.config.doScaleVariance: 

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

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

311 backgrounds = afwMath.BackgroundList() 

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

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

314 sources = detections.sources 

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

316 for bg in detections.background: 

317 backgrounds.append(bg) 

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

319 

320 def _cropToExactBinning(self, exposure, patchInfo): 

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

322 background binning. 

323 

324 Parameters 

325 ---------- 

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

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

328 patchInfo : `lsst.skymap.PatchInfo` 

329 Description of the patch geometry. 

330 

331 Returns 

332 ------- 

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

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

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

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

337 contained by the patch outer bounding box. 

338 

339 Raises 

340 ------ 

341 ValueError 

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

343 of the background bin size. 

344 """ 

345 bbox = patchInfo.getInnerBBox() 

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

347 raise ValueError( 

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

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

350 ) 

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

352 raise ValueError( 

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

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

355 ) 

356 outer_bbox = patchInfo.getOuterBBox() 

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

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

359 bbox.grow( 

360 Extent2I( 

361 n_bins_grow_x*self.detection.background.binSizeX, 

362 n_bins_grow_y*self.detection.background.binSizeY, 

363 ) 

364 ) 

365 assert outer_bbox.contains(bbox) 

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

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

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

369 return exposure[bbox] 

370 

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

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

373 

374 Parameters 

375 ---------- 

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

377 Exposure that the background should correspond to. 

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

379 Description of the patch geometry. Only needed if 

380 `~DetectCoaddSourceConfig.forceExactBinning` is `True`. 

381 

382 Returns 

383 ------- 

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

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

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

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

388 subtraction. 

389 """ 

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

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

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

393 # "stats image" with NaNs and NO_DATA. 

394 if self.config.forceExactBinning: 

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

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

397 exposure.mask[:, :] = 0 

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

399 for bg, *_ in background: 

400 stats = bg.getStatsImage() 

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

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

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

404 return background 

405 

406 

407class MeasureMergedCoaddSourcesConnections( 

408 PipelineTaskConnections, 

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

410 defaultTemplates={ 

411 "inputCoaddName": "deep", 

412 "outputCoaddName": "deep", 

413 "deblendedCatalog": "deblendedFlux", 

414 }, 

415 deprecatedTemplates={ 

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

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

418 }, 

419): 

420 inputSchema = cT.InitInput( 

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

422 name="{inputCoaddName}Coadd_deblendedFlux_schema", 

423 storageClass="SourceCatalog" 

424 ) 

425 outputSchema = cT.InitOutput( 

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

427 name="{inputCoaddName}Coadd_meas_schema", 

428 storageClass="SourceCatalog" 

429 ) 

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

431 refCat = cT.PrerequisiteInput( 

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

433 name="ref_cat", 

434 storageClass="SimpleCatalog", 

435 dimensions=("skypix",), 

436 deferLoad=True, 

437 multiple=True, 

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

439 ) 

440 exposure = cT.Input( 

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

442 name="{inputCoaddName}Coadd_calexp", 

443 storageClass="ExposureF", 

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

445 ) 

446 exposure_cells = cT.Input( 

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

448 name="{inputCoaddName}CoaddCell", 

449 storageClass="MultipleCellCoadd", 

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

451 ) 

452 background = cT.Input( 

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

454 name="{inputCoaddName}Coadd_calexp_background", 

455 storageClass="Background", 

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

457 ) 

458 skyMap = cT.Input( 

459 doc="SkyMap to use in processing", 

460 name=BaseSkyMap.SKYMAP_DATASET_TYPE_NAME, 

461 storageClass="SkyMap", 

462 dimensions=("skymap",), 

463 ) 

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

465 visitCatalogs = cT.Input( 

466 doc="Deprecated and unused.", 

467 name="src", 

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

469 storageClass="SourceCatalog", 

470 multiple=True, 

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

472 ) 

473 sourceTableHandles = cT.Input( 

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

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

476 "PSF flags."), 

477 name="sourceTable_visit", 

478 storageClass="ArrowAstropy", 

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

480 multiple=True, 

481 deferLoad=True, 

482 ) 

483 finalizedSourceTableHandles = cT.Input( 

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

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

486 name="finalized_src_table", 

487 storageClass="ArrowAstropy", 

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

489 multiple=True, 

490 deferLoad=True, 

491 ) 

492 finalVisitSummaryHandles = cT.Input( 

493 doc="Final visit summary table", 

494 name="finalVisitSummary", 

495 storageClass="ExposureCatalog", 

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

497 multiple=True, 

498 deferLoad=True, 

499 ) 

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

501 inputCatalog = cT.Input( 

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

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

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

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

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

507 "be 'mergeDet'"), 

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

509 storageClass="SourceCatalog", 

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

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

512 ) 

513 scarletCatalog = cT.Input( 

514 doc="Catalogs produced by multiband deblending", 

515 name="{inputCoaddName}Coadd_deblendedCatalog", 

516 storageClass="SourceCatalog", 

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

518 ) 

519 scarletModels = cT.Input( 

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

521 name="{inputCoaddName}Coadd_scarletModelData", 

522 storageClass="LsstScarletModelData", 

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

524 ) 

525 outputSources = cT.Output( 

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

527 name="{outputCoaddName}Coadd_meas", 

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

529 storageClass="SourceCatalog", 

530 ) 

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

532 matchResult = cT.Output( 

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

534 name="{outputCoaddName}Coadd_measMatch", 

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

536 storageClass="Catalog", 

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

538 ) 

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

540 denormMatches = cT.Output( 

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

542 "doWriteMatchesDenormalized", 

543 name="{outputCoaddName}Coadd_measMatchFull", 

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

545 storageClass="Catalog", 

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

547 ) 

548 

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

550 super().__init__(config=config) 

551 del self.visitCatalogs 

552 if not config.doPropagateFlags: 

553 del self.sourceTableHandles 

554 del self.finalizedSourceTableHandles 

555 else: 

556 # Check for types of flags required. 

557 if not config.propagateFlags.source_flags: 

558 del self.sourceTableHandles 

559 if not config.propagateFlags.finalized_source_flags: 

560 del self.finalizedSourceTableHandles 

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

562 if config.inputCatalog == "deblendedCatalog": 

563 del self.inputCatalog 

564 if not config.doAddFootprints: 

565 del self.scarletModels 

566 else: 

567 del self.deblendedCatalog 

568 del self.scarletModels 

569 

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

571 if not config.doMatchSources: 

572 del self.refCat 

573 del self.matchResult 

574 

575 if not config.doWriteMatchesDenormalized: 

576 del self.denormMatches 

577 

578 if config.useCellCoadds: 

579 del self.exposure 

580 else: 

581 del self.exposure_cells 

582 del self.background 

583 

584 

585class MeasureMergedCoaddSourcesConfig(PipelineTaskConfig, 

586 pipelineConnections=MeasureMergedCoaddSourcesConnections): 

587 """Configuration parameters for the MeasureMergedCoaddSourcesTask 

588 """ 

589 inputCatalog = ChoiceField( 

590 dtype=str, 

591 default="deblendedCatalog", 

592 allowed={ 

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

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

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

596 }, 

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

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

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

600 ) 

601 doAddFootprints = Field(dtype=bool, 

602 default=True, 

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

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

605 "otherwise this should be False.") 

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

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

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

609 "model footprints.") 

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

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

612 "saving to disk. " 

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

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

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

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

617 doPropagateFlags = Field( 

618 dtype=bool, default=True, 

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

620 ) 

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

622 doMatchSources = Field( 

623 dtype=bool, 

624 default=False, 

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

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

627 ) 

628 match = ConfigurableField( 

629 target=DirectMatchTask, 

630 doc="Matching to reference catalog", 

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

632 ) 

633 doWriteMatchesDenormalized = Field( 

634 dtype=bool, 

635 default=False, 

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

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

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

639 ) 

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

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

642 checkUnitsParseStrict = Field( 

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

644 dtype=str, 

645 default="raise", 

646 ) 

647 doApCorr = Field( 

648 dtype=bool, 

649 default=True, 

650 doc="Apply aperture corrections" 

651 ) 

652 applyApCorr = ConfigurableField( 

653 target=ApplyApCorrTask, 

654 doc="Subtask to apply aperture corrections" 

655 ) 

656 doRunCatalogCalculation = Field( 

657 dtype=bool, 

658 default=True, 

659 doc='Run catalogCalculation task' 

660 ) 

661 catalogCalculation = ConfigurableField( 

662 target=CatalogCalculationTask, 

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

664 ) 

665 

666 hasFakes = Field( 

667 dtype=bool, 

668 default=False, 

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

670 ) 

671 idGenerator = SkyMapIdGeneratorConfig.make_field() 

672 

673 @property 

674 def refObjLoader(self): 

675 return self.match.refObjLoader 

676 

677 def setDefaults(self): 

678 super().setDefaults() 

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

680 'base_Variance', 

681 'base_LocalPhotoCalib', 

682 'base_LocalWcs'] 

683 

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

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

686 # still set flags with it here. 

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

688 'INEXACT_PSF'] 

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

690 'INEXACT_PSF'] 

691 

692 def validate(self): 

693 super().validate() 

694 

695 if not self.doMatchSources and self.doWriteMatchesDenormalized: 

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

697 

698 

699class MeasureMergedCoaddSourcesTask(PipelineTask): 

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

701 

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

703 measurement in each coadd. 

704 

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

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

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

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

709 

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

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

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

713 propagated to the coadd sources. 

714 

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

716 catalog. 

717 

718 After MeasureMergedCoaddSourcesTask has been run on multiple coadds, we 

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

720 processing procedure will merge these measurements into a suitable catalog 

721 for driving forced photometry. 

722 

723 Parameters 

724 ---------- 

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

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

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

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

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

730 An instance of ReferenceObjectLoader that supplies an external reference 

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

732 requiring a reference catalog are disabled. 

733 initInputs : `dict`, optional 

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

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

736 **kwargs 

737 Additional keyword arguments. 

738 """ 

739 

740 _DefaultName = "measureCoaddSources" 

741 ConfigClass = MeasureMergedCoaddSourcesConfig 

742 

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

744 **kwargs): 

745 super().__init__(**kwargs) 

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

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

748 if initInputs is not None: 

749 schema = initInputs['inputSchema'].schema 

750 if schema is None: 

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

752 self.schemaMapper = afwTable.SchemaMapper(schema) 

753 self.schemaMapper.addMinimalSchema(schema) 

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

755 self.algMetadata = PropertyList() 

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

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

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

759 if self.config.doMatchSources: 

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

761 if self.config.doPropagateFlags: 

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

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

764 if self.config.doApCorr: 

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

766 if self.config.doRunCatalogCalculation: 

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

768 

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

770 

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

772 inputs = butlerQC.get(inputRefs) 

773 

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

775 if self.config.doMatchSources: 

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

777 inputs.pop('refCat'), 

778 name=self.config.connections.refCat, 

779 config=self.config.refObjLoader, 

780 log=self.log) 

781 self.match.setRefObjLoader(refObjLoader) 

782 

783 if self.config.useCellCoadds: 

784 multiple_cell_coadd = inputs.pop("exposure_cells") 

785 stitched_coadd = multiple_cell_coadd.stitch() 

786 exposure = stitched_coadd.asExposure() 

787 background = inputs.pop("background") 

788 exposure.image -= background.getImage() 

789 

790 ccdInputs = stitched_coadd.ccds 

791 apCorrMap = stitched_coadd.ap_corr_map 

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

793 else: 

794 exposure = inputs.pop("exposure") 

795 # Set psfcache 

796 # move this to run after gen2 deprecation 

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

798 

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

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

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

802 

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

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

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

806 

807 # Transform inputCatalog 

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

809 sources = afwTable.SourceCatalog(table) 

810 # Load the correct input catalog 

811 if "scarletCatalog" in inputs: 

812 inputCatalog = inputs.pop("scarletCatalog") 

813 catalogRef = inputRefs.scarletCatalog 

814 else: 

815 inputCatalog = inputs.pop("inputCatalog") 

816 catalogRef = inputRefs.inputCatalog 

817 sources.extend(inputCatalog, self.schemaMapper) 

818 del inputCatalog 

819 # Add the HeavyFootprints to the deblended sources 

820 if self.config.doAddFootprints: 

821 modelData = inputs.pop('scarletModels') 

822 if self.config.doConserveFlux: 

823 imageForRedistribution = exposure 

824 else: 

825 imageForRedistribution = None 

826 updateCatalogFootprints( 

827 modelData=modelData, 

828 catalog=sources, 

829 band=band, 

830 imageForRedistribution=imageForRedistribution, 

831 removeScarletData=True, 

832 updateFluxColumns=True, 

833 ) 

834 table = sources.getTable() 

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

836 

837 skyMap = inputs.pop('skyMap') 

838 tractNumber = catalogRef.dataId['tract'] 

839 tractInfo = skyMap[tractNumber] 

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

841 skyInfo = Struct( 

842 skyMap=skyMap, 

843 tractInfo=tractInfo, 

844 patchInfo=patchInfo, 

845 wcs=tractInfo.getWcs(), 

846 bbox=patchInfo.getOuterBBox() 

847 ) 

848 

849 if self.config.doPropagateFlags: 

850 if "sourceTableHandles" in inputs: 

851 sourceTableHandles = inputs.pop("sourceTableHandles") 

852 sourceTableHandleDict = {handle.dataId["visit"]: handle for handle in sourceTableHandles} 

853 else: 

854 sourceTableHandleDict = None 

855 if "finalizedSourceTableHandles" in inputs: 

856 finalizedSourceTableHandles = inputs.pop("finalizedSourceTableHandles") 

857 finalizedSourceTableHandleDict = {handle.dataId["visit"]: handle 

858 for handle in finalizedSourceTableHandles} 

859 else: 

860 finalizedSourceTableHandleDict = None