Coverage for python / lsst / pipe / tasks / matchBackgrounds.py: 0%
213 statements
« prev ^ index » next coverage.py v7.13.5, created at 2026-04-25 08:39 +0000
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/>.
22__all__ = [
23 "MatchBackgroundsConnections",
24 "MatchBackgroundsConfig",
25 "MatchBackgroundsTask",
26 "ChooseReferenceVisitConfig",
27 "ChooseReferenceVisitTask",
28]
30import numpy as np
32from lsst.afw.image import ImageF, MaskedImageF
33from lsst.afw.math import (
34 MEANSQUARE,
35 NPOINT,
36 VARIANCE,
37 ApproximateControl,
38 BackgroundControl,
39 BackgroundMI,
40 StatisticsControl,
41 makeBackground,
42 makeStatistics,
43 stringToInterpStyle,
44 stringToStatisticsProperty,
45 stringToUndersampleStyle,
46)
47from lsst.pex.config import ChoiceField, Config, ConfigField, ConfigurableField, Field, RangeField
48from lsst.pipe.base import PipelineTask, PipelineTaskConfig, PipelineTaskConnections, Struct, Task
49from lsst.pipe.base.connectionTypes import Input, Output
50from lsst.pipe.tasks.tractBackground import TractBackground, TractBackgroundConfig
51from lsst.skymap import BaseSkyMap
52from lsst.utils.timer import timeMethod
55class ChooseReferenceVisitConfig(Config):
57 tractBgModel = ConfigField(
58 dtype=TractBackgroundConfig,
59 doc="Background model for the entire tract",
60 )
61 bestRefWeightChi2 = RangeField(
62 dtype=float,
63 doc="Mean background goodness of fit statistic weight when calculating the best reference exposure. "
64 "Higher weights prefer exposures with flatter backgrounds. Ignored when ref visit supplied.",
65 default=0.3,
66 min=0.0,
67 max=1.0,
68 )
69 bestRefWeightVariance = RangeField(
70 dtype=float,
71 doc="Image variance weight when calculating the best reference exposure. "
72 "Higher weights prefers exposures with low image variances. Ignored when ref visit supplied.",
73 default=0.3,
74 min=0.0,
75 max=1.0,
76 )
77 bestRefWeightGlobalCoverage = RangeField(
78 dtype=float,
79 doc="Global coverage weight (total number of valid pixels) when calculating the best reference "
80 "exposure. Higher weights prefer exposures with high coverage. Ignored when ref visit supplied.",
81 default=0.4,
82 min=0.0,
83 max=1.0,
84 )
85 gridStatistic = ChoiceField(
86 dtype=str,
87 doc="Type of statistic to estimate pixel value for the grid points",
88 default="MEANCLIP",
89 allowed={"MEAN": "mean", "MEDIAN": "median", "MEANCLIP": "clipped mean"},
90 )
91 undersampleStyle = ChoiceField(
92 dtype=str,
93 doc="Behavior if there are too few points in the grid for requested interpolation style",
94 default="REDUCE_INTERP_ORDER",
95 allowed={
96 "THROW_EXCEPTION": "throw an exception if there are too few points.",
97 "REDUCE_INTERP_ORDER": "use an interpolation style with a lower order.",
98 },
99 )
100 interpStyle = ChoiceField(
101 dtype=str,
102 doc="Algorithm to interpolate the background values; ignored if ``usePolynomial=True``."
103 "Maps to an enum; see afw.math.Background for more information.",
104 default="AKIMA_SPLINE",
105 allowed={
106 "CONSTANT": "Use a single constant value.",
107 "LINEAR": "Use linear interpolation.",
108 "NATURAL_SPLINE": "A cubic spline with zero second derivative at endpoints.",
109 "AKIMA_SPLINE": "A higher-level non-linear spline that is more robust to outliers.",
110 "NONE": "No background estimation is to be attempted.",
111 },
112 )
113 numSigmaClip = Field[int](
114 doc="Sigma for outlier rejection. Ignored if ``gridStatistic != 'MEANCLIP'``.",
115 default=3,
116 )
117 numIter = Field[int](
118 doc="Number of iterations of outlier rejection. Ignored if ``gridStatistic != 'MEANCLIP'``.",
119 default=3,
120 )
123class ChooseReferenceVisitTask(Task):
124 """Select a reference visit from a list of visits by minimizing a cost
125 function
127 Notes
128 -----
129 The reference exposure is chosen from the list of science exposures by
130 minimizing a cost function that penalizes high background complexity
131 (divergence from a plane), high variance, and low global coverage.
132 The cost function is a weighted sum of these three metrics.
133 The weights are set by the config parameters:
135 - ``bestRefWeightChi2``
136 - ``bestRefWeightVariance``
137 - ``bestRefWeightGlobalCoverage``
138 """
140 # TODO: this is only used to select reference visits, so may become
141 # redundant if we find a means to build reference images instead.
142 ConfigClass = ChooseReferenceVisitConfig
143 config: ChooseReferenceVisitConfig
145 def __init__(self, *args, **kwargs):
146 super().__init__(**kwargs)
147 # Fits on binned images only; masking controlled in tractBackground.py
148 self.statsFlag = stringToStatisticsProperty(self.config.gridStatistic)
149 self.statsCtrl = StatisticsControl()
150 self.statsCtrl.setNanSafe(True)
151 self.statsCtrl.setNumSigmaClip(self.config.numSigmaClip)
152 self.statsCtrl.setNumIter(self.config.numIter)
153 self.stringToInterpStyle = stringToInterpStyle(self.config.interpStyle)
154 self.undersampleStyle = stringToUndersampleStyle(self.config.undersampleStyle)
156 @timeMethod
157 def _makeTractBackgrounds(self, warps, skyMap):
158 """Create full tract models of all visit backgrounds
160 Parameters
161 ----------
162 warps : `list`[`~lsst.daf.butler.DeferredDatasetHandle`]
163 List of warped exposures (of type `~lsst.afw.image.ExposureF`).
164 This is ordered by patch ID, then by visit ID
165 skyMap : `lsst.skyMap.SkyMap`
166 SkyMap for deriving the patch/tract dimensions
168 Returns
169 -------
170 visitTractBackrounds : `dict` [`int`, `TractBackground`]
171 Models of full tract backgrounds for all visits, in nanojanskies.
172 Accessed by visit ID.
173 """
174 # First, separate warps by visit
175 visits = np.unique([i.dataId["visit"] for i in warps])
177 # Then build background models for each visit and store
178 visitTractBackgrounds = {}
179 for i in range(len(visits)):
180 visitWarpDDFs = [j for j in warps if j.dataId["visit"] == visits[i]]
181 # Set up empty full tract background model object
182 bgModelBase = TractBackground(
183 config=self.config.tractBgModel, skymap=skyMap, tract=warps[0].dataId["tract"]
184 )
186 bgModels = []
187 for warp in visitWarpDDFs:
188 msg = "Constructing FFP background model for visit %d using %d patches"
189 self.log.debug(
190 msg,
191 visits[i],
192 len(visitWarpDDFs),
193 )
194 workingWarp = warp.get()
196 bgModel = bgModelBase.clone()
197 bgModel.addWarp(workingWarp)
198 bgModels.append(bgModel)
200 # Merge warp models to make a single full tract background model
201 for bgModel, warp in zip(bgModels, visitWarpDDFs):
202 msg = (
203 "Patch %d: Merging %d unmasked pixels (%.1f%s of detector area) into full tract "
204 "background model"
205 )
206 self.log.debug(
207 msg,
208 warp.dataId["patch"],
209 bgModel._numbers.getArray().sum(),
210 100 * bgModel._numbers.getArray().sum() / workingWarp.getBBox().getArea(),
211 "%",
212 )
213 bgModelBase += bgModel
215 visitTractBackgrounds[visits[i]] = bgModelBase
217 return visitTractBackgrounds
219 @timeMethod
220 def _defineWarps(self, visitTractBackgrounds):
221 """Define the reference visit
223 This method calculates an appropriate reference visit from the
224 supplied full tract visit backgrounds by minimizing a cost function
225 that penalizes high background complexity (divergence from a plane),
226 high variance within the warps comprising the visit, and low global
227 coverage.
229 Parameters
230 ----------
231 visitTractBackrounds : `dict` [`int`, `TractBackground`]
232 Models of full tract backgrounds for all visits, in nanojanskies.
233 Accessed by visit ID.
234 Returns
235 -------
236 refVisId : `int`
237 ID of the reference visit.
238 """
239 # Extract mean/var/npoints for each visit background model
240 fitChi2s = [] # Background goodness of fit
241 visitVars = [] # Variance of original image
242 fitNPointsGlobal = [] # Global coverage
243 visits = [] # To ensure dictionary key order is correct
244 for vis in visitTractBackgrounds:
245 visits.append(vis)
246 # Fit a polynomial model to the full tract plane
247 tractBg = visitTractBackgrounds[vis].getStatsImage()
248 tractVar = visitTractBackgrounds[vis].getVarianceImage()
249 fitBg, _ = fitBackground(tractBg, self.statsCtrl, self.statsFlag, self.config.undersampleStyle)
250 # Weight the approximation fit by the binned image variance
251 fitBg.getStatsImage().variance = tractVar
253 # Return an approximation to the background
254 approxCtrl = ApproximateControl(ApproximateControl.CHEBYSHEV, 1, 1, True)
255 fitApprox = fitBg.getApproximate(approxCtrl, self.undersampleStyle)
257 fitBgSub = MaskedImageF(ImageF(tractBg.array - fitApprox.getImage().array))
258 bad_mask_bit_mask = fitBgSub.mask.getPlaneBitMask("BAD")
259 fitBgSub.mask.array[np.isnan(fitBgSub.image.array)] = bad_mask_bit_mask
261 fitStats = makeStatistics(fitBgSub.image, fitBgSub.mask, VARIANCE | NPOINT, self.statsCtrl)
263 good = (fitBgSub.mask.array.astype(int) & bad_mask_bit_mask) == 0
264 dof = len(good[good]) - 6 # Assuming eq. of plane
265 fitChi2 = (
266 np.nansum(fitBgSub.image.array[good] ** 2 / fitBg.getStatsImage().variance.array[good]) / dof
267 )
269 visitVar = np.nanmean(fitBg.getStatsImage().variance.array[good])
270 fitNPointGlobal, _ = fitStats.getResult(NPOINT)
271 fitChi2s.append(fitChi2)
272 visitVars.append(visitVar)
273 fitNPointsGlobal.append(int(fitNPointGlobal))
275 self.log.info(
276 "Sci exp. visit %d; BG fit Chi^2=%.2e, var=%.2f nJy, nPoints global=%d",
277 vis,
278 fitChi2,
279 visitVar,
280 fitNPointGlobal,
281 )
282 # Normalize mean/var/npoints to range from 0 to 1
283 fitChi2sFrac = np.array(fitChi2s) / np.nanmax(fitChi2s)
284 fitVarsFrac = np.array(visitVars) / np.nanmax(visitVars)
285 fitNPointsGlobalFrac = np.nanmin(fitNPointsGlobal) / np.array(fitNPointsGlobal)
287 # Calculate cost function values
288 costFunctionVals = self.config.bestRefWeightChi2 * fitChi2sFrac
289 costFunctionVals += self.config.bestRefWeightVariance * fitVarsFrac
290 costFunctionVals += self.config.bestRefWeightGlobalCoverage * fitNPointsGlobalFrac
292 ind = np.nanargmin(costFunctionVals)
293 refVisitId = visits[ind]
294 self.log.info("Using best reference visit %d", refVisitId)
296 return refVisitId
299class MatchBackgroundsConnections(
300 PipelineTaskConnections,
301 dimensions=("skymap", "tract", "band"),
302 defaultTemplates={
303 "warpType": "psf_matched",
304 "warpTypeSuffix": "",
305 },
306):
307 # TODO: add connection for warped backgroundToPhotometricRatio
308 warps = Input(
309 doc=("Warps used to construct a list of matched backgrounds."),
310 name="{warpType}_warp",
311 storageClass="ExposureF",
312 dimensions=("skymap", "tract", "patch", "visit"),
313 deferLoad=True,
314 multiple=True,
315 )
316 skyMap = Input(
317 doc="Input definition of geometry/bbox and projection/wcs for warped exposures",
318 name=BaseSkyMap.SKYMAP_DATASET_TYPE_NAME,
319 storageClass="SkyMap",
320 dimensions=("skymap",),
321 )
322 backgroundInfoList = Output(
323 doc="List of differential backgrounds, with goodness of fit params.",
324 name="{warpType}_warp_background_diff", # TODO: settle on appropriate name
325 dimensions=("skymap", "tract", "visit", "patch"),
326 storageClass="Background",
327 multiple=True,
328 )
329 matchedImageList = Output(
330 doc="List of background-matched warps.",
331 name="{warpType}_warp_background_matched", # TODO: settle on appropriate name
332 storageClass="ExposureF",
333 dimensions=("skymap", "tract", "visit", "patch"),
334 multiple=True,
335 )
338class MatchBackgroundsConfig(PipelineTaskConfig, pipelineConnections=MatchBackgroundsConnections):
340 refWarpVisit = Field[int](
341 doc="Reference visit ID. If None, the best visit is chosen using the list of visits.",
342 optional=True,
343 )
344 tractBgModel = ConfigField(
345 dtype=TractBackgroundConfig,
346 doc="Background model for the entire tract",
347 )
348 gridStatistic = ChoiceField(
349 dtype=str,
350 doc="Type of statistic to estimate pixel value for the grid points",
351 default="MEANCLIP",
352 allowed={"MEAN": "mean", "MEDIAN": "median", "MEANCLIP": "clipped mean"},
353 )
354 undersampleStyle = ChoiceField(
355 dtype=str,
356 doc="Behavior if there are too few points in the grid for requested interpolation style",
357 default="REDUCE_INTERP_ORDER",
358 allowed={
359 "THROW_EXCEPTION": "throw an exception if there are too few points.",
360 "REDUCE_INTERP_ORDER": "use an interpolation style with a lower order.",
361 },
362 )
363 interpStyle = ChoiceField(
364 dtype=str,
365 doc="Algorithm to interpolate the background values. "
366 "Maps to an enum; see afw.math.Background for more information.",
367 default="AKIMA_SPLINE",
368 allowed={
369 "CONSTANT": "Use a single constant value.",
370 "LINEAR": "Use linear interpolation.",
371 "NATURAL_SPLINE": "A cubic spline with zero second derivative at endpoints.",
372 "AKIMA_SPLINE": "A higher-level non-linear spline that is more robust to outliers.",
373 "NONE": "No background estimation is to be attempted.",
374 },
375 )
376 numSigmaClip = Field[int](
377 doc="Sigma for outlier rejection. Ignored if ``gridStatistic != 'MEANCLIP'``.",
378 default=3,
379 )
380 numIter = Field[int](
381 doc="Number of iterations of outlier rejection. Ignored if ``gridStatistic != 'MEANCLIP'``.",
382 default=3,
383 )
384 reference = ConfigurableField(
385 target=ChooseReferenceVisitTask, doc="Choose reference visit to match backgrounds to"
386 )
388 def validate(self):
389 if self.undersampleStyle == "INCREASE_NXNYSAMPLE":
390 raise RuntimeError("Invalid undersampleStyle: Polynomial fitting not implemented.")
393class MatchBackgroundsTask(PipelineTask):
394 """Match the backgrounds of a list of warped exposures to a reference
396 Attributes
397 ----------
398 config : `MatchBackgroundsConfig`
399 Configuration for this task.
400 statsCtrl : `~lsst.afw.math.StatisticsControl`
401 Statistics control object.
403 Notes
404 -----
405 This task is a part of the background subtraction pipeline. It matches the
406 backgrounds of a list of warped science exposures to that of a reference
407 image.
408 """
410 ConfigClass = MatchBackgroundsConfig
411 config: MatchBackgroundsConfig
412 _DefaultName = "matchBackgrounds"
414 def __init__(self, *args, **kwargs):
415 super().__init__(**kwargs)
416 # Fits on binned images only; masking controlled in tractBackground.py
417 self.statsFlag = stringToStatisticsProperty(self.config.gridStatistic)
418 self.statsCtrl = StatisticsControl()
419 self.statsCtrl.setNanSafe(True)
420 self.statsCtrl.setNumSigmaClip(self.config.numSigmaClip)
421 self.statsCtrl.setNumIter(self.config.numIter)
422 self.stringToInterpStyle = stringToInterpStyle(self.config.interpStyle)
423 self.undersampleStyle = stringToUndersampleStyle(self.config.undersampleStyle)
425 self.makeSubtask("reference")
427 @timeMethod
428 def run(self, warps, skyMap):
429 """Match the backgrounds of a list of warped exposures to the same
430 patches in a reference visit
432 A reference visit ID will be chosen automatically if none is supplied.
434 Parameters
435 ----------
436 warps : `list`[`~lsst.afw.image.Exposure`]
437 List of warped science exposures to be background-matched.
438 skyMap : `lsst.skyMap.SkyMap`
439 SkyMap for deriving the patch/tract dimensions.
441 Returns
442 -------
443 result : `~lsst.afw.math.BackgroundList`, `~lsst.afw.image.Exposure`
444 Differential background models and associated background-matched
445 images.
446 """
447 # TODO: include warped backgroundToPhotometricRatio correction
448 if (numExp := len(warps)) < 1:
449 self.log.warning("No exposures found! Returning empty lists.")
450 return Struct(backgroundInfoList=[], matchedImageList=[])
452 if self.config.refWarpVisit is None:
453 # Build FFP BG models of each visit
454 visitTractBgs = self.reference._makeTractBackgrounds(warps, skyMap)
455 # Choose a reference visit using those
456 refVisId = self.reference._defineWarps(visitTractBgs)
457 else:
458 self.log.info("Using user-supplied reference visit %d", self.config.refWarpVisit)
459 refVisId = self.config.refWarpVisit
461 self.log.info("Matching %d Exposures", numExp)
463 backgroundInfoList, matchedImageList = self.matchBackgrounds(warps, skyMap, refVisId)
465 return Struct(backgroundInfoList=backgroundInfoList, matchedImageList=matchedImageList)
467 @timeMethod
468 def _makeTractDifferenceBackgrounds(self, warps, skyMap, refVisitId):
469 """Create full tract models of all difference image backgrounds
470 (reference visit - visit)
472 Parameters
473 ----------
474 warps : `list`[`~lsst.daf.butler.DeferredDatasetHandle`]
475 List of warped exposures (of type `~lsst.afw.image.ExposureF`).
476 This is ordered by patch ID, then by visit ID
477 skyMap : `lsst.skyMap.SkyMap`
478 SkyMap for deriving the patch/tract dimensions.
479 refVisitId : `int`
480 Visit ID number for the chosen reference visit.
482 Returns
483 -------
484 visitTractDifferenceBackrounds : `dict` [`int`, `TractBackground`]
485 Models of full tract backgrounds for all difference images
486 (reference visit - visit), in nanojanskies.
487 Accessed by visit ID.
488 """
489 # First, separate warps by visit
490 visits = np.unique([i.dataId["visit"] for i in warps])
492 # Then build difference image background models for each visit & store
493 visitTractDifferenceBackgrounds = {}
494 for i in range(len(visits)):
495 visitWarpDDFs = [j for j in warps if j.dataId["visit"] == visits[i]]
496 refWarpDDFs = [j for j in warps if j.dataId["visit"] == refVisitId]
497 refPatches = [j.dataId["patch"] for j in refWarpDDFs]
498 # Set up empty full tract background model object
499 bgModelBase = TractBackground(
500 config=self.config.tractBgModel, skymap=skyMap, tract=warps[0].dataId["tract"]
501 )
503 bgModels = []
504 for warp in visitWarpDDFs:
505 msg = "Constructing FFP background model for reference visit %d - visit %d using %d patches"
506 self.log.debug(
507 msg,
508 refVisitId,
509 visits[i],
510 len(visitWarpDDFs),
511 )
512 workingWarp = warp.get()
514 patchId = warp.dataId["patch"]
515 # On no overlap between working warp and reference visit, set
516 # the image to all NaN
517 try:
518 idx = refPatches.index(patchId)
519 refWarp = refWarpDDFs[idx].get()
520 except ValueError:
521 refWarp = workingWarp.clone()
522 refWarp.image += np.nan
523 workingWarp.image.array = refWarp.image.array - workingWarp.image.array
525 bgModel = bgModelBase.clone()
526 bgModel.addWarp(workingWarp)
527 bgModels.append(bgModel)
529 # Merge warp difference models to make a single full tract
530 # background difference model
531 for bgModel, warp in zip(bgModels, visitWarpDDFs):
532 msg = (
533 "Patch %d: Merging %d unmasked pixels (%.1f%s of detector area) into full tract "
534 "difference background model"
535 )
536 self.log.debug(
537 msg,
538 warp.dataId["patch"],
539 bgModel._numbers.getArray().sum(),
540 100 * bgModel._numbers.getArray().sum() / workingWarp.getBBox().getArea(),
541 "%",
542 )
543 bgModelBase += bgModel
545 # Fit full tract background to generate offset image
546 if visits[i] != refVisitId:
547 bgModelImage = bgModelBase.getStatsImage()
548 # Note: this just extrapolates into regions of no overlap
549 # between reference and visit
550 bkgd, _ = fitBackground(
551 bgModelImage, self.statsCtrl, self.statsFlag, self.config.undersampleStyle
552 )
553 try:
554 bkgdImage = bkgd.getImageF(self.config.interpStyle, self.config.undersampleStyle)
555 except Exception as e:
556 e.add_note(f"on image {warp.dataId}")
557 raise
558 # Calculate RMS and MSE of fit and print as log
559 resids = ImageF(bgModelImage.array - bkgdImage.array)
560 rms = np.sqrt(np.nanmean(resids.array**2))
561 mse = makeStatistics(resids, MEANSQUARE, self.statsCtrl).getValue()
563 self.log.info(
564 "Visit %d; difference BG fit RMS=%.2f nJy, matched MSE=%.2f nJy",
565 visits[i],
566 rms,
567 mse,
568 )
569 # Replace binned difference image w/best-fit model.
570 # Resetting numbers to override interpolation
571 bgModelBase._numbers.array[:] = 1e6 # Arbitrarily large value
572 bgModelBase._values.array = bkgdImage.array * bgModelBase._numbers.array
574 visitTractDifferenceBackgrounds[visits[i]] = bgModelBase
576 return visitTractDifferenceBackgrounds
578 @timeMethod
579 def matchBackgrounds(self, warps, skyMap, refVisitId):
580 """Match science visit's background level to that of reference visit
582 Process creates binned images of the full focal plane (in tract
583 coordinates) for all visit IDs, subtracts each from a similarly
584 binned FFP reference image, then generates TractBackground
585 objects.
587 The TractBackground objects representing the difference image
588 backgrounds are then used to generate 'offset' images for each warp
589 comprising the full science exposure visit, which are then added to
590 each warp to match the background to that of the reference visit at the
591 warp's location within the tract.
593 Best practice uses `psf_matched_warp` images without the
594 detections mask plane set. When using `direct_warp`, sources may bias
595 the difference image background estimation. Mask planes are set in
596 TractBackgroundConfig.
598 Fit diagnostics are also calculated and returned.
600 Parameters
601 ----------
602 warps : `list`[`~lsst.daf.butler.DeferredDatasetHandle`]
603 List of warped exposures (of type `~lsst.afw.image.ExposureF`).
604 This is ordered by patch ID, then by visit ID
605 skyMap : `lsst.skyMap.SkyMap`
606 SkyMap for deriving the patch/tract dimensions.
607 refVisitId : `int`
608 Chosen reference visit ID to match to.
610 Returns
611 -------
612 backgroundInfoList : `list`[`TractBackground`]
613 List of all difference image backgrounds used to match to reference
614 visit warps, in nanojanskies.
615 matchedImageList : `list`[`~lsst.afw.image.ExposureF`]
616 List of all background-matched warps, in nanojanskies.
617 """
618 visits = np.unique([i.dataId["visit"] for i in warps])
619 self.log.info("Processing %d visits", len(visits))
621 backgroundInfoList = []
622 matchedImageList = []
623 diffTractBackgrounds = self._makeTractDifferenceBackgrounds(warps, skyMap, refVisitId)
625 # Reference visit doesn't need an offset image, so use all 0's
626 im = warps[0].get() # Use arbitrary image as base
627 bkgd = diffTractBackgrounds[refVisitId].toWarpBackground(im)
628 blank = bkgd.getImage()
629 blank *= 0
631 for warp in warps:
632 visId = warp.dataId["visit"]
633 if visId == refVisitId:
634 backgroundInfoList.append(bkgd) # Just append a 0 image
635 matchedImageList.append(warp.get())
636 continue
637 self.log.info(
638 "Matching background of %s to same patch in visit %s",
639 warp.dataId,
640 refVisitId,
641 )
642 im = warp.get()
643 maskIm = im.getMaskedImage()
644 tractBg = diffTractBackgrounds[visId]
645 diffModel = tractBg.toWarpBackground(im)
646 bkgdIm = diffModel.getImage()
647 maskIm.image += bkgdIm
649 backgroundInfoList.append(diffModel)
650 matchedImageList.append(im)
652 return backgroundInfoList, matchedImageList
655def fitBackground(
656 warp: MaskedImageF, statsCtrl, statsFlag, undersampleStyle
657) -> tuple[BackgroundMI, BackgroundControl]:
658 """Generate a simple binned background masked image for warped or other
659 data
661 Parameters
662 ----------
663 warp: `~lsst.afw.image.MaskedImageF`
664 Warped exposure for which to estimate background.
665 statsCtrl : `~lsst.afw.math.StatisticsControl`
666 Statistics control object.
667 statsFlag : `~lsst.afw.math.Property`
668 Statistics control type.
669 undersampleStyle : `str`
670 Behavior if there are too few points in the grid for requested
671 interpolation style.
673 Returns
674 -------
675 bkgd: `~lsst.afw.math.BackgroundMI`
676 Background model of masked warp.
677 bgCtrl: `~lsst.afw.math.BackgroundControl`
678 Background control object.
679 """
680 # This only accesses pre-binned images, so no scaling necessary
681 nx = warp.getWidth()
682 ny = warp.getHeight()
684 bgCtrl = BackgroundControl(nx, ny, statsCtrl, statsFlag)
685 bgCtrl.setUndersampleStyle(undersampleStyle)
686 bkgd = makeBackground(warp, bgCtrl)
688 return bkgd, bgCtrl