Coverage for python / lsst / pipe / tasks / repair.py: 21%
112 statements
« prev ^ index » next coverage.py v7.13.5, created at 2026-05-06 08:52 +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__ = ["RepairConfig", "RepairTask", "TooManyCosmicRays"]
24import lsst.pex.config as pexConfig
25import lsst.afw.math as afwMath
26import lsst.afw.detection as afwDet
27import lsst.meas.algorithms as measAlg
28import lsst.pipe.base as pipeBase
29from lsstDebug import getDebugFrame
30import lsst.afw.display as afwDisplay
31from lsst.pipe.tasks.interpImage import InterpImageTask
32from lsst.utils.timer import timeMethod
33from lsst.pex.exceptions import LengthError
36class TooManyCosmicRays(pipeBase.AlgorithmError):
37 """Raised if the cosmic ray task fails with too many cosmics.
39 Parameters
40 ----------
41 maxCosmicRays : `int`
42 Maximum number of cosmic rays allowed.
43 """
44 def __init__(self, maxCosmicRays, **kwargs):
45 msg = f"Cosmic ray task found more than {maxCosmicRays} cosmics."
46 self.msg = msg
47 self._metadata = kwargs
48 super().__init__(msg, **kwargs)
49 self._metadata["maxCosmicRays"] = maxCosmicRays
51 def __str__(self):
52 # Exception doesn't handle **kwargs, so we need a custom str.
53 return f"{self.msg}: {self.metadata}"
55 @property
56 def metadata(self):
57 for key, value in self._metadata.items():
58 if not (isinstance(value, int) or isinstance(value, float) or isinstance(value, str)):
59 raise TypeError(f"{key} is of type {type(value)}, but only (int, float, str) are allowed.")
60 return self._metadata
63class RepairConfig(pexConfig.Config):
64 doInterpolate = pexConfig.Field(
65 dtype=bool,
66 doc="Interpolate over defects? (ignored unless you provide a list of defects)",
67 default=True,
68 )
69 doCosmicRay = pexConfig.Field(
70 dtype=bool,
71 doc="Find and mask out cosmic rays?",
72 default=True,
73 )
74 cosmicray = pexConfig.ConfigField(
75 dtype=measAlg.FindCosmicRaysConfig,
76 doc="Options for finding and masking cosmic rays",
77 )
78 interp = pexConfig.ConfigurableField(
79 target=InterpImageTask,
80 doc="Interpolate over bad image pixels",
81 )
83 def setDefaults(self):
84 self.interp.useFallbackValueAtEdge = True
85 self.interp.fallbackValueType = "MEANCLIP"
86 self.interp.negativeFallbackAllowed = True
89class RepairTask(pipeBase.Task):
90 """Repair an exposures defects and cosmic rays via interpolation.
92 This task operates on an lsst.afw.image.Exposure in place to
93 interpolate over a set of `~lsst.meas.algorithms.Defect` objects.
95 It will also, optionally, find and interpolate any cosmic rays in the lsst.afw.image.Exposure.
97 Notes
98 -----
99 Debugging:
100 The available debug variables in RepairTask are:
102 display :
103 A dictionary containing debug point names as keys with frame number as value. Valid keys are:
104 repair.before :
105 display image before any repair is done
106 repair.after :
107 display image after cosmic ray and defect correction
108 displayCR :
109 If True, display the exposure on ds9's frame 1 and overlay bounding boxes around detects CRs.
111 To investigate the pipe_tasks_repair_Debug, put something like
113 .. code-block :: none
115 import lsstDebug
116 def DebugInfo(name):
117 di = lsstDebug.getInfo(name) # N.b. lsstDebug.Info(name) would call us recursively
118 if name == "lsst.pipe.tasks.repair":
119 di.display = {'repair.before':2, 'repair.after':3}
120 di.displayCR = True
121 return di
123 lsstDebug.Info = DebugInfo
124 into your debug.py file and run runRepair.py with the --debug flag.
126 Conversion notes:
127 Display code should be updated once we settle on a standard way of controlling what is displayed.
128 """
130 ConfigClass = RepairConfig
131 _DefaultName = "repair"
133 def __init__(self, **kwargs):
134 pipeBase.Task.__init__(self, **kwargs)
135 if self.config.doInterpolate:
136 self.makeSubtask("interp")
138 @timeMethod
139 def run(self, exposure, defects=None, keepCRs=None):
140 """Repair an Exposure's defects and cosmic rays.
142 Parameters
143 ----------
144 exposure : `lsst.afw.image.Exposure`
145 Exposure must have a valid Psf.
146 Modified in place.
147 defects : `lsst.meas.algorithms.DefectListT` or `None`, optional
148 If `None`, do no defect correction.
149 keepCRs : `Unknown` or `None`, optional
150 Don't interpolate over the CR pixels (defer to ``RepairConfig`` if `None`).
152 Raises
153 ------
154 AssertionError
155 Raised if any of the following occur:
156 - No exposure provided.
157 - The object provided as exposure evaluates to False.
158 - No PSF provided.
159 - The Exposure has no associated Psf.
160 """
161 assert exposure, "No exposure provided"
162 psf = exposure.getPsf()
163 assert psf, "No PSF provided"
165 frame = getDebugFrame(self._display, "repair.before")
166 if frame:
167 afwDisplay.Display(frame).mtv(exposure)
169 if defects is not None and self.config.doInterpolate:
170 self.interp.run(exposure, defects=defects)
172 if self.config.doCosmicRay:
173 self.cosmicRay(exposure, keepCRs=keepCRs)
175 frame = getDebugFrame(self._display, "repair.after")
176 if frame:
177 afwDisplay.Display(frame).mtv(exposure)
179 def cosmicRay(self, exposure, keepCRs=None):
180 """Mask cosmic rays.
182 Parameters
183 ----------
184 exposure : `lsst.afw.image.Exposure`
185 Exposure to process.
186 keepCRs : `Unknown` or `None`, optional
187 Don't interpolate over the CR pixels (defer to ``pex_config`` if `None`).
188 """
189 import lsstDebug
190 display = lsstDebug.Info(__name__).display
191 displayCR = lsstDebug.Info(__name__).displayCR
193 assert exposure, "No exposure provided"
194 psf = exposure.getPsf()
195 assert psf, "No psf provided"
197 # Blow away old mask
198 try:
199 mask = exposure.getMaskedImage().getMask()
200 crBit = mask.getMaskPlane("CR")
201 mask.clearMaskPlane(crBit)
202 except Exception:
203 pass
205 exposure0 = exposure # initial value of exposure
206 binSize = self.config.cosmicray.background.binSize
207 nx, ny = exposure.getWidth()/binSize, exposure.getHeight()/binSize
208 # Treat constant background as a special case to avoid the extra complexity in calling
209 # measAlg.SubtractBackgroundTask().
210 if nx*ny <= 1:
211 medianBg = afwMath.makeStatistics(exposure.getMaskedImage(), afwMath.MEDIAN).getValue()
212 modelBg = None
213 else:
214 # make a deep copy of the exposure before subtracting its background,
215 # because this routine is only allowed to modify the exposure by setting mask planes
216 # and interpolating over defects, not changing the background level
217 exposure = exposure.Factory(exposure, True)
218 subtractBackgroundTask = measAlg.SubtractBackgroundTask(config=self.config.cosmicray.background)
219 modelBg = subtractBackgroundTask.run(exposure).background
220 medianBg = 0.0
222 if keepCRs is None:
223 keepCRs = self.config.cosmicray.keepCRs
224 try:
225 try:
226 crs = measAlg.findCosmicRays(
227 exposure.getMaskedImage(),
228 psf,
229 medianBg,
230 pexConfig.makePropertySet(self.config.cosmicray),
231 keepCRs,
232 )
233 except LengthError:
234 raise TooManyCosmicRays(self.config.cosmicray.nCrPixelMax) from None
236 if modelBg:
237 # Add back background image
238 img = exposure.getMaskedImage()
239 img += modelBg.getImageF()
240 del img
241 # Replace original image with CR subtracted image
242 exposure0.setMaskedImage(exposure.getMaskedImage())
244 except Exception:
245 if display:
246 afwDisplay.Display().mtv(exposure0, title="Failed CR")
247 raise
249 num = 0
250 if crs is not None:
251 mask = exposure0.getMaskedImage().getMask()
252 crBit = mask.getPlaneBitMask("CR")
253 afwDet.setMaskFromFootprintList(mask, crs, crBit)
254 num = len(crs)
256 if display and displayCR:
257 disp = afwDisplay.Display()
258 disp.incrDefaultFrame()
259 disp.mtv(exposure0, title="Post-CR")
261 with disp.Buffering():
262 for cr in crs:
263 afwDisplay.utils.drawBBox(cr.getBBox(), borderWidth=0.55)
265 text = "kept" if keepCRs else "interpolated over"
266 self.log.info("Identified and %s %s cosmic rays.", text, num)
267 self.metadata["cosmic_ray_count"] = num