Coverage for python / lsst / pipe / tasks / skyCorrection.py: 17%

220 statements  

« prev     ^ index     » next       coverage.py v7.13.1, created at 2026-01-09 01:51 -0800

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 

22from __future__ import annotations 

23 

24__all__ = ["SkyCorrectionTask", "SkyCorrectionConfig"] 

25 

26import warnings 

27 

28import numpy as np 

29 

30from lsst.afw.image import ExposureF, makeMaskedImage 

31from lsst.afw.math import BackgroundMI, binImage 

32from lsst.daf.butler import DeferredDatasetHandle 

33from lsst.pex.config import Config, ConfigField, ConfigurableField, Field 

34from lsst.pipe.base import PipelineTask, PipelineTaskConfig, PipelineTaskConnections, Struct 

35from lsst.pipe.base.connectionTypes import Input, Output, PrerequisiteInput 

36from lsst.pipe.tasks.background import ( 

37 FocalPlaneBackground, 

38 FocalPlaneBackgroundConfig, 

39 MaskObjectsTask, 

40 SkyMeasurementTask, 

41) 

42from lsst.pipe.tasks.visualizeVisit import VisualizeMosaicExpConfig, VisualizeMosaicExpTask 

43 

44 

45def _reorderAndPadList(inputList, inputKeys, outputKeys, padWith=None): 

46 """Match the order of one list to another, padding if necessary. 

47 

48 Parameters 

49 ---------- 

50 inputList : `list` 

51 List to be reordered and padded. Elements can be any type. 

52 inputKeys : iterable 

53 Iterable of values to be compared with outputKeys. 

54 Length must match `inputList`. 

55 outputKeys : iterable 

56 Iterable of values to be compared with inputKeys. 

57 padWith : 

58 Any value to be inserted where one of inputKeys is not in outputKeys. 

59 

60 Returns 

61 ------- 

62 outputList : `list` 

63 Copy of inputList reordered per outputKeys and padded with `padWith` 

64 so that the length matches length of outputKeys. 

65 """ 

66 outputList = [] 

67 for outputKey in outputKeys: 

68 if outputKey in inputKeys: 

69 outputList.append(inputList[inputKeys.index(outputKey)]) 

70 else: 

71 outputList.append(padWith) 

72 return outputList 

73 

74 

75class SkyCorrectionConnections(PipelineTaskConnections, dimensions=("instrument", "visit")): 

76 calExps = Input( 

77 doc="Background-subtracted calibrated exposures.", 

78 name="calexp", 

79 multiple=True, 

80 deferLoad=True, 

81 storageClass="ExposureF", 

82 dimensions=["instrument", "visit", "detector"], 

83 ) 

84 calBkgs = Input( 

85 doc="Subtracted backgrounds for input calibrated exposures.", 

86 name="calexpBackground", 

87 multiple=True, 

88 storageClass="Background", 

89 dimensions=["instrument", "visit", "detector"], 

90 ) 

91 backgroundToPhotometricRatioHandles = Input( 

92 doc="Ratio of a background-flattened image to a photometric-flattened image. " 

93 "Only used if doApplyFlatBackgroundRatio is True.", 

94 multiple=True, 

95 name="background_to_photometric_ratio", 

96 storageClass="Image", 

97 dimensions=["instrument", "visit", "detector"], 

98 deferLoad=True, 

99 ) 

100 skyFrames = PrerequisiteInput( 

101 doc="Calibration sky frames.", 

102 name="sky", 

103 multiple=True, 

104 deferLoad=True, 

105 storageClass="ExposureF", 

106 dimensions=["instrument", "physical_filter", "detector"], 

107 isCalibration=True, 

108 ) 

109 camera = PrerequisiteInput( 

110 doc="Input camera.", 

111 name="camera", 

112 storageClass="Camera", 

113 dimensions=["instrument"], 

114 isCalibration=True, 

115 ) 

116 skyCorr = Output( 

117 doc="Sky correction data, to be subtracted from the calibrated exposures.", 

118 name="skyCorr", 

119 multiple=True, 

120 storageClass="Background", 

121 dimensions=["instrument", "visit", "detector"], 

122 ) 

123 calExpMosaic = Output( 

124 doc="Full focal plane mosaicked image of the sky corrected calibrated exposures.", 

125 name="calexp_skyCorr_visit_mosaic", 

126 storageClass="ImageF", 

127 dimensions=["instrument", "visit"], 

128 ) 

129 calBkgMosaic = Output( 

130 doc="Full focal plane mosaicked image of the sky corrected calibrated exposure backgrounds.", 

131 name="calexpBackground_skyCorr_visit_mosaic", 

132 storageClass="ImageF", 

133 dimensions=["instrument", "visit"], 

134 ) 

135 

136 def __init__(self, *, config: "SkyCorrectionConfig | None" = None): 

137 super().__init__(config=config) 

138 assert config is not None 

139 if not config.doSky: 

140 del self.skyFrames 

141 if not config.doApplyFlatBackgroundRatio: 

142 del self.backgroundToPhotometricRatioHandles 

143 

144 

145class SkyCorrectionConfig(PipelineTaskConfig, pipelineConnections=SkyCorrectionConnections): 

146 doApplyFlatBackgroundRatio = Field( 

147 dtype=bool, 

148 default=False, 

149 doc="This should be True if the input image was processed with an illumination correction.", 

150 ) 

151 maskObjects = ConfigurableField( 

152 target=MaskObjectsTask, 

153 doc="Mask Objects", 

154 ) 

155 doMaskObjects = Field( 

156 dtype=bool, 

157 default=True, 

158 doc="Iteratively mask objects to find good sky?", 

159 ) 

160 bgModel1 = ConfigField( 

161 dtype=FocalPlaneBackgroundConfig, 

162 doc="Initial background model, prior to sky frame subtraction", 

163 ) 

164 undoBgModel1 = Field( 

165 dtype=bool, 

166 default=False, 

167 doc="If True, adds back initial background model after sky and removes bgModel1 from the list", 

168 ) 

169 sky = ConfigurableField( 

170 target=SkyMeasurementTask, 

171 doc="Sky measurement", 

172 ) 

173 doSky = Field( 

174 dtype=bool, 

175 default=True, 

176 doc="Do sky frame subtraction?", 

177 ) 

178 bgModel2 = ConfigField( 

179 dtype=FocalPlaneBackgroundConfig, 

180 doc="Final (cleanup) background model, after sky frame subtraction", 

181 ) 

182 doBgModel2 = Field( 

183 dtype=bool, 

184 default=True, 

185 doc="Do final (cleanup) background model subtraction, after sky frame subtraction?", 

186 ) 

187 binning = Field( 

188 dtype=int, 

189 default=8, 

190 doc="Binning factor for constructing full focal plane '*_camera' output datasets", 

191 ) 

192 

193 def setDefaults(self): 

194 Config.setDefaults(self) 

195 self.bgModel2.doSmooth = True 

196 self.bgModel2.minFrac = 0.3 

197 self.bgModel2.xSize = 256 

198 self.bgModel2.ySize = 256 

199 self.bgModel2.smoothScale = 1.0 

200 

201 def validate(self): 

202 super().validate() 

203 if self.undoBgModel1 and not self.doSky and not self.doBgModel2: 

204 raise ValueError("If undoBgModel1 is True, task requires at least one of doSky or doBgModel2.") 

205 

206 

207class SkyCorrectionTask(PipelineTask): 

208 """Perform a full focal plane sky correction.""" 

209 

210 ConfigClass = SkyCorrectionConfig 

211 _DefaultName = "skyCorr" 

212 

213 def __init__(self, *args, **kwargs): 

214 super().__init__(**kwargs) 

215 self.makeSubtask("sky") 

216 self.makeSubtask("maskObjects") 

217 

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

219 # Sort input/output connections by detector ID, padding where 

220 # necessary, to ensure that all detectors are processed consistently. 

221 # Detector IDs are defined from the intersection of calExps, calBkgs, 

222 # and optionally skyFrames and backgroundToPhotometricRatioHandles. 

223 # This resolves potential missing data issues when processing a visit 

224 # that contains only partial inputs. 

225 calExpOrder = {ref.dataId["detector"] for ref in inputRefs.calExps} 

226 calBkgOrder = {ref.dataId["detector"] for ref in inputRefs.calBkgs} 

227 detectorOrder = calExpOrder & calBkgOrder 

228 if self.config.doApplyFlatBackgroundRatio: 

229 ratioOrder = {ref.dataId["detector"] for ref in inputRefs.backgroundToPhotometricRatioHandles} 

230 detectorOrder &= ratioOrder 

231 if self.config.doSky: 

232 skyFrameOrder = {ref.dataId["detector"] for ref in inputRefs.skyFrames} 

233 detectorOrder &= skyFrameOrder 

234 detectorOrder = sorted(detectorOrder) 

235 inputRefs.calExps = _reorderAndPadList( 

236 inputRefs.calExps, [ref.dataId["detector"] for ref in inputRefs.calExps], detectorOrder 

237 ) 

238 inputRefs.calBkgs = _reorderAndPadList( 

239 inputRefs.calBkgs, [ref.dataId["detector"] for ref in inputRefs.calBkgs], detectorOrder 

240 ) 

241 # Only attempt to fetch sky frames if they are going to be applied. 

242 if self.config.doSky: 

243 inputRefs.skyFrames = _reorderAndPadList( 

244 inputRefs.skyFrames, [ref.dataId["detector"] for ref in inputRefs.skyFrames], detectorOrder 

245 ) 

246 else: 

247 inputRefs.skyFrames = [] 

248 # Only attempt to fetch flat ratios if they are going to be applied. 

249 if self.config.doApplyFlatBackgroundRatio: 

250 inputRefs.backgroundToPhotometricRatioHandles = _reorderAndPadList( 

251 inputRefs.backgroundToPhotometricRatioHandles, 

252 [ref.dataId["detector"] for ref in inputRefs.backgroundToPhotometricRatioHandles], 

253 detectorOrder, 

254 ) 

255 else: 

256 inputRefs.backgroundToPhotometricRatioHandles = [] 

257 outputRefs.skyCorr = _reorderAndPadList( 

258 outputRefs.skyCorr, [ref.dataId["detector"] for ref in outputRefs.skyCorr], detectorOrder 

259 ) 

260 inputs = butlerQC.get(inputRefs) 

261 outputs = self.run(**inputs) 

262 butlerQC.put(outputs, outputRefs) 

263 

264 def run(self, calExps, calBkgs, skyFrames, camera, backgroundToPhotometricRatioHandles=[]): 

265 """Perform sky correction on a visit. 

266 

267 The original visit-level background is first restored to the calibrated 

268 exposure and the existing background model is inverted in-place. If 

269 doMaskObjects is True, the mask map associated with this exposure will 

270 be iteratively updated (over nIter loops) by re-estimating the 

271 background each iteration and redetecting footprints. 

272 

273 An initial full focal plane sky subtraction (bgModel1) will take place 

274 prior to scaling and subtracting the sky frame. 

275 

276 If doSky is True, the sky frame will be scaled to the flux in the input 

277 visit. 

278 

279 If doBgModel2 is True, a final full focal plane sky subtraction will 

280 take place after the sky frame has been subtracted. 

281 

282 The first N elements of the returned skyCorr will consist of inverted 

283 elements of the calexpBackground model (i.e., subtractive). All 

284 subsequent elements appended to skyCorr thereafter will be additive 

285 such that, when skyCorr is subtracted from a calexp, the net result 

286 will be to undo the initial per-detector background solution and then 

287 apply the skyCorr model thereafter. Adding skyCorr to a 

288 calexpBackground will effectively negate the calexpBackground, 

289 returning only the additive background components of the skyCorr 

290 background model. 

291 

292 Parameters 

293 ---------- 

294 calExps : `list` [`lsst.afw.image.ExposureF`] 

295 Detector calibrated exposure images for the visit. 

296 calBkgs : `list` [`lsst.afw.math.BackgroundList`] 

297 Detector background lists matching the calibrated exposures. 

298 skyFrames : `list` [`lsst.afw.image.ExposureF`] 

299 Sky frame calibration data for the input detectors. 

300 camera : `lsst.afw.cameraGeom.Camera` 

301 Camera matching the input data to process. 

302 backgroundToPhotometricRatioHandles : 

303 `list` [`lsst.daf.butler.DeferredDatasetHandle`], optional 

304 Deferred dataset handles pointing to the Background to photometric 

305 ratio images for the input detectors. 

306 

307 Returns 

308 ------- 

309 results : `Struct` containing: 

310 skyFrameScale : `float` 

311 Scale factor applied to the sky frame. 

312 skyCorr : `list` [`lsst.afw.math.BackgroundList`] 

313 Detector-level sky correction background lists. 

314 calExpMosaic : `lsst.afw.image.ExposureF` 

315 Visit-level mosaic of the sky corrected data, binned. 

316 Analogous to `calexp - skyCorr`. 

317 calBkgMosaic : `lsst.afw.image.ExposureF` 

318 Visit-level mosaic of the sky correction background, binned. 

319 Analogous to `calexpBackground + skyCorr`. 

320 """ 

321 # Process each detector separately and merge into bgModel1. 

322 # This avoids storing every full-res image in-memory at once. 

323 bgModel1 = FocalPlaneBackground.fromCamera(self.config.bgModel1, camera) 

324 detectors = [] 

325 masks = [] 

326 skyCorrs = [] 

327 bgModel1Indices = [] 

328 if not self.config.doApplyFlatBackgroundRatio: 

329 backgroundToPhotometricRatioHandles = [None] * len(calExps) 

330 for calExpHandle, calBkg, backgroundToPhotometricRatioHandle in zip( 

331 calExps, calBkgs, backgroundToPhotometricRatioHandles 

332 ): 

333 calExp = self._getCalExp( 

334 calExpHandle, backgroundToPhotometricRatioHandle=backgroundToPhotometricRatioHandle 

335 ) 

336 detectors.append(calExp.getDetector()) 

337 

338 # Restore original background in-place; optionally refine mask maps 

339 _ = self._restoreOriginalBackgroundRefineMask(calExp, calBkg) 

340 masks.append(calExp.mask) 

341 skyCorrs.append(calBkg) # Contains only the inverted original background elements at this stage 

342 bgModel1Indices.append(len(calBkg)) # Index of the original background element 

343 

344 # Make a background model for the image, using bgModel1 configs 

345 bgModel1Detector = FocalPlaneBackground.fromCamera(self.config.bgModel1, camera) 

346 bgModel1Detector.addCcd(calExp) 

347 bgModel1.merge(bgModel1Detector) 

348 self.log.info( 

349 "Detector %d: Merged %d unmasked pixels (%.1f%s of detector area) into initial BG model", 

350 calExp.getDetector().getId(), 

351 bgModel1Detector._numbers.getArray().sum(), 

352 100 * bgModel1Detector._numbers.getArray().sum() / calExp.getBBox().getArea(), 

353 "%", 

354 ) 

355 

356 # Validate bgModel1 

357 self._validateBgModel("bgModel1", bgModel1, self.config.bgModel1) 

358 

359 # Update skyCorrs with new bgModel1 background elements 

360 for detector, skyCorr in zip(detectors, skyCorrs): 

361 with warnings.catch_warnings(): 

362 warnings.filterwarnings("ignore", "invalid value encountered") 

363 calBkgElement = bgModel1.toCcdBackground(detector, detector.getBBox()) 

364 skyCorr.append(calBkgElement[0]) 

365 

366 # Fit a scaled sky frame to all input exposures 

367 skyFrameScale = None 

368 if self.config.doSky: 

369 skyFrameScale = self._fitSkyFrame( 

370 calExps, masks, skyCorrs, skyFrames, backgroundToPhotometricRatioHandles 

371 ) 

372 

373 # Remove the initial background model (bgModel1) from every skyCorr 

374 if self.config.undoBgModel1: 

375 for skyCorr, bgModel1Index in zip(skyCorrs, bgModel1Indices): 

376 skyCorr._backgrounds.pop(bgModel1Index) 

377 self.log.info( 

378 "Initial background models (bgModel1s) have been removed from all skyCorr background lists", 

379 ) 

380 

381 # Bin exposures, generate full-fp bg, map to CCDs and subtract in-place 

382 if self.config.doBgModel2: 

383 bgModel2 = FocalPlaneBackground.fromCamera(self.config.bgModel2, camera) 

384 for calExpHandle, mask, skyCorr, backgroundToPhotometricRatioHandle in zip( 

385 calExps, masks, skyCorrs, backgroundToPhotometricRatioHandles 

386 ): 

387 calExp = self._getCalExp(calExpHandle, mask, skyCorr, backgroundToPhotometricRatioHandle) 

388 

389 # Make a background model for the image, using bgModel2 configs 

390 bgModel2Detector = FocalPlaneBackground.fromCamera(self.config.bgModel2, camera) 

391 bgModel2Detector.addCcd(calExp) 

392 bgModel2.merge(bgModel2Detector) 

393 self.log.info( 

394 "Detector %d: Merged %d unmasked pixels (%.1f%s of detector area) into final BG model", 

395 calExp.getDetector().getId(), 

396 bgModel2Detector._numbers.getArray().sum(), 

397 100 * bgModel2Detector._numbers.getArray().sum() / calExp.getBBox().getArea(), 

398 "%", 

399 ) 

400 

401 # Validate bgModel2 

402 self._validateBgModel("bgModel2", bgModel2, self.config.bgModel2) 

403 

404 # Update skyCorrs with new bgModel2 background elements 

405 for detector, skyCorr in zip(detectors, skyCorrs): 

406 with warnings.catch_warnings(): 

407 warnings.filterwarnings("ignore", "invalid value encountered") 

408 calBkgElement = bgModel2.toCcdBackground(detector, detector.getBBox()) 

409 skyCorr.append(calBkgElement[0]) 

410 

411 # Make camera-level mosaics of bg subtracted calexps and subtracted bgs 

412 calExpsBinned = [] 

413 calBkgsBinned = [] 

414 for calExpHandle, mask, skyCorr, backgroundToPhotometricRatioHandle, bgModel1Index in zip( 

415 calExps, masks, skyCorrs, backgroundToPhotometricRatioHandles, bgModel1Indices 

416 ): 

417 calExp = self._getCalExp(calExpHandle, mask, skyCorr, backgroundToPhotometricRatioHandle) 

418 

419 skyCorrExtra = skyCorr.clone() # new skyCorr elements created in this task 

420 skyCorrExtra._backgrounds = skyCorrExtra._backgrounds[bgModel1Index:] 

421 skyCorrExtraMI = makeMaskedImage(skyCorrExtra.getImage()) 

422 skyCorrExtraMI.setMask(calExp.getMask()) 

423 

424 calExpsBinned.append(binImage(calExp.getMaskedImage(), self.config.binning)) 

425 calBkgsBinned.append(binImage(skyCorrExtraMI, self.config.binning)) 

426 

427 mosConfig = VisualizeMosaicExpConfig() 

428 mosConfig.binning = self.config.binning 

429 mosTask = VisualizeMosaicExpTask(config=mosConfig) 

430 detectorIds = [detector.getId() for detector in detectors] 

431 calExpMosaic = mosTask.run(calExpsBinned, camera, inputIds=detectorIds).outputData 

432 calBkgMosaic = mosTask.run(calBkgsBinned, camera, inputIds=detectorIds).outputData 

433 

434 return Struct( 

435 skyFrameScale=skyFrameScale, 

436 skyCorr=skyCorrs, 

437 calExpMosaic=calExpMosaic, 

438 calBkgMosaic=calBkgMosaic, 

439 ) 

440 

441 def _getCalExp(self, calExpHandle, mask=None, skyCorr=None, backgroundToPhotometricRatioHandle=None): 

442 """Get a calexp from a DeferredDatasetHandle, and optionally apply an 

443 updated mask and skyCorr. 

444 

445 Parameters 

446 ---------- 

447 calExpHandle : `~lsst.afw.image.ExposureF` 

448 | `lsst.daf.butler.DeferredDatasetHandle` 

449 Either the image exposure data or a handle to the calexp dataset. 

450 mask : `lsst.afw.image.Mask`, optional 

451 Mask to apply to the calexp. 

452 skyCorr : `lsst.afw.math.BackgroundList`, optional 

453 Background list to subtract from the calexp. 

454 

455 Returns 

456 ------- 

457 calExp : `lsst.afw.image.ExposureF` 

458 The calexp with the mask and skyCorr applied. 

459 """ 

460 if isinstance(calExpHandle, DeferredDatasetHandle): 

461 calExp: ExposureF = calExpHandle.get() 

462 else: 

463 # Here we clone the imaging data to avoid modifying data which is 

464 # used in downstream processing. 

465 calExp: ExposureF = calExpHandle.clone() 

466 

467 # Convert from background-flattened to photometric-flattened images 

468 # Note: remember to convert back to background-flattened images 

469 # if science images are to be output by this task. 

470 if self.config.doApplyFlatBackgroundRatio: 

471 if not backgroundToPhotometricRatioHandle: 

472 raise ValueError( 

473 "A list of backgroundToPhotometricRatioHandles must be supplied if " 

474 "config.doApplyFlatBackgroundRatio=True.", 

475 ) 

476 ratioImage = backgroundToPhotometricRatioHandle.get() 

477 calExp.maskedImage *= ratioImage 

478 self.log.info( 

479 "Detector %d: Converted background-flattened image to a photometric-flattened image", 

480 calExp.getDetector().getId(), 

481 ) 

482 

483 if mask is not None: 

484 calExp.setMask(mask) 

485 if skyCorr is not None: 

486 image = calExp.getMaskedImage() 

487 image -= skyCorr.getImage() 

488 

489 return calExp 

490 

491 def _getSkyFrame(self, skyFrameHandle): 

492 """Get a calexp from a DeferredDatasetHandle, and optionally apply an 

493 updated mask and skyCorr. 

494 

495 Parameters 

496 ---------- 

497 skyFrameHandle : `lsst.daf.butler.DeferredDatasetHandle` 

498 Either the sky frame data or a handle to the sky frame dataset. 

499 

500 Returns 

501 ------- 

502 skyFrame : `lsst.afw.image.ExposureF` 

503 The calexp with the mask and skyCorr applied. 

504 """ 

505 if isinstance(skyFrameHandle, DeferredDatasetHandle): 

506 skyFrame: ExposureF = skyFrameHandle.get() 

507 else: 

508 skyFrame: ExposureF = skyFrameHandle 

509 return skyFrame 

510 

511 def _restoreOriginalBackgroundRefineMask(self, calExp, calBkg): 

512 """Restore original background to a calexp and invert the related 

513 background model; optionally refine the mask plane. 

514 

515 The original visit-level background is restored to the calibrated 

516 exposure and the existing background model is inverted in-place. If 

517 doMaskObjects is True, the mask map associated with the exposure will 

518 be iteratively updated (over nIter loops) by re-estimating the 

519 background each iteration and redetecting footprints. 

520 

521 The background model modified in-place in this method will comprise the 

522 first N elements of the skyCorr dataset type, i.e., these N elements 

523 are the inverse of the calexpBackground model. All subsequent elements 

524 appended to skyCorr will be additive such that, when skyCorr is 

525 subtracted from a calexp, the net result will be to undo the initial 

526 per-detector background solution and then apply the skyCorr model 

527 thereafter. Adding skyCorr to a calexpBackground will effectively 

528 negate the calexpBackground, returning only the additive background 

529 components of the skyCorr background model. 

530 

531 Parameters 

532 ---------- 

533 calExp : `lsst.afw.image.ExposureF` 

534 Detector level calexp image. 

535 calBkg : `lsst.afw.math.BackgroundList` 

536 Detector level background lists associated with the calexp. 

537 

538 Returns 

539 ------- 

540 calExp : `lsst.afw.image.ExposureF` 

541 The calexp with the originally subtracted background restored. 

542 skyCorrBase : `lsst.afw.math.BackgroundList` 

543 The inverted original background models; the genesis for skyCorr. 

544 """ 

545 image = calExp.getMaskedImage() 

546 

547 # Invert all elements of the existing bg model; restore in calexp 

548 for calBkgElement in calBkg: 

549 statsImage = calBkgElement[0].getStatsImage() 

550 statsImage *= -1 

551 skyCorrBase = calBkg.getImage() 

552 image -= skyCorrBase 

553 

554 stats = np.nanpercentile(skyCorrBase.array, [50, 75, 25]) 

555 self.log.info( 

556 "Detector %d: Original background restored (BG median = %.1f counts, BG IQR = %.1f counts)", 

557 calExp.getDetector().getId(), 

558 -stats[0], 

559 np.subtract(*stats[1:]), 

560 ) 

561 

562 # Iteratively subtract bg, re-detect sources, and add bg back on 

563 if self.config.doMaskObjects: 

564 maskFrac0 = 1 - np.sum(calExp.mask.array == 0) / calExp.mask.array.size 

565 self.maskObjects.findObjects(calExp) 

566 maskFrac1 = 1 - np.sum(calExp.mask.array == 0) / calExp.mask.array.size 

567 

568 self.log.info( 

569 "Detector %d: Iterative source detection and mask growth has increased masked area by %.1f%%", 

570 calExp.getDetector().getId(), 

571 (100 * (maskFrac1 - maskFrac0)), 

572 ) 

573 

574 return calExp, skyCorrBase 

575 

576 def _validateBgModel(self, bgModelID, bgModel, config): 

577 """Check that the background model contains enough valid superpixels, 

578 and raise a useful error if not. 

579 

580 Parameters 

581 ---------- 

582 bgModelID : `str` 

583 Identifier for the background model. 

584 bgModel : `~lsst.pipe.tasks.background.FocalPlaneBackground` 

585 Background model to check. 

586 config : `~lsst.pipe.tasks.background.FocalPlaneBackgroundConfig` 

587 Configuration used to create the background model. 

588 """ 

589 bgModelArray = bgModel._numbers.getArray() 

590 spArea = (config.xSize / config.pixelSize) * (config.ySize / config.pixelSize) 

591 self.log.info( 

592 "%s: FP background model constructed using %.2f x %.2f mm (%d x %d pixel) superpixels", 

593 bgModelID, 

594 config.xSize, 

595 config.ySize, 

596 int(config.xSize / config.pixelSize), 

597 int(config.ySize / config.pixelSize), 

598 ) 

599 self.log.info( 

600 "%s: Pixel data exists in %d of %d superpixels; the most populated superpixel is %.1f%% filled", 

601 bgModelID, 

602 np.sum(bgModelArray > 0), 

603 bgModelArray.size, 

604 100 * np.max(bgModelArray) / spArea, 

605 ) 

606 

607 thresh = config.minFrac * spArea 

608 if np.all(bgModelArray < thresh): 

609 raise RuntimeError( 

610 f"No background model superpixels are more than {100*config.minFrac}% filled. " 

611 "Try decreasing the minFrac configuration parameter, optimizing the subset of detectors " 

612 "being processed, or increasing the number of detectors being processed." 

613 ) 

614 

615 with warnings.catch_warnings(): 

616 warnings.filterwarnings("ignore", r"invalid value encountered") 

617 stats = np.nanpercentile(bgModel.getStatsImage().array, [50, 75, 25]) 

618 self.log.info( 

619 "%s: FP BG median = %.1f counts, FP BG IQR = %.1f counts", 

620 bgModelID, 

621 stats[0], 

622 np.subtract(*stats[1:]), 

623 ) 

624 

625 def _fitSkyFrame(self, calExps, masks, skyCorrs, skyFrames, backgroundToPhotometricRatioHandles): 

626 """Determine the full focal plane sky frame scale factor relative to 

627 an input list of calibrated exposures. 

628 

629 This method measures the sky frame scale on all inputs, resulting in 

630 values equal to the background method solveScales(). 

631 

632 Input skyCorrs are updated in-place. 

633 

634 Parameters 

635 ---------- 

636 calExps : `list` [`lsst.afw.image.ExposureF`] 

637 Calibrated exposures to be background subtracted. 

638 masks : `list` [`lsst.afw.image.Mask`] 

639 Masks associated with the input calibrated exposures. 

640 skyCorrs : `list` [`lsst.afw.math.BackgroundList`] 

641 Background lists associated with the input calibrated exposures. 

642 skyFrames : `list` [`lsst.afw.image.ExposureF`] 

643 Sky frame calibration data for the input detectors. 

644 

645 Returns 

646 ------- 

647 scale : `float` 

648 Fitted scale factor applied to the sky frame. 

649 """ 

650 skyBkgs = [] 

651 scales = [] 

652 for calExpHandle, mask, skyCorr, skyFrameHandle, backgroundToPhotometricRatioHandle in zip( 

653 calExps, masks, skyCorrs, skyFrames, backgroundToPhotometricRatioHandles 

654 ): 

655 calExp = self._getCalExp(calExpHandle, mask, skyCorr, backgroundToPhotometricRatioHandle) 

656 skyFrame = self._getSkyFrame(skyFrameHandle) 

657 skyBkg = self.sky.exposureToBackground(skyFrame) 

658 del skyFrame # Free up memory 

659 skyBkgs.append(skyBkg) 

660 # Return a tuple of gridded image and sky frame clipped means 

661 samples = self.sky.measureScale(calExp.getMaskedImage(), skyBkg) 

662 scales.append(samples) 

663 scale = self.sky.solveScales(scales) 

664 for skyCorr, skyBkg in zip(skyCorrs, skyBkgs): 

665 bgData = list(skyBkg[0]) 

666 bg = bgData[0] 

667 statsImage = bg.getStatsImage().clone() 

668 statsImage *= scale 

669 newBg = BackgroundMI(bg.getImageBBox(), statsImage) 

670 newBgData = [newBg] + bgData[1:] 

671 skyCorr.append(newBgData) 

672 self.log.info("Sky frame subtracted with a scale factor of %.5f", scale) 

673 return scale