LSST Applications g04e9c324dd+8c5ae1fdc5,g134cb467dc+1b3060144d,g18429d2f64+f642bf4753,g199a45376c+0ba108daf9,g1fd858c14a+2dcf163641,g262e1987ae+7b8c96d2ca,g29ae962dfc+3bd6ecb08a,g2cef7863aa+aef1011c0b,g35bb328faa+8c5ae1fdc5,g3fd5ace14f+53e1a9e7c5,g4595892280+fef73a337f,g47891489e3+2efcf17695,g4d44eb3520+642b70b07e,g53246c7159+8c5ae1fdc5,g67b6fd64d1+2efcf17695,g67fd3c3899+b70e05ef52,g74acd417e5+317eb4c7d4,g786e29fd12+668abc6043,g87389fa792+8856018cbb,g89139ef638+2efcf17695,g8d7436a09f+3be3c13596,g8ea07a8fe4+9f5ccc88ac,g90f42f885a+a4e7b16d9b,g97be763408+ad77d7208f,g9dd6db0277+b70e05ef52,ga681d05dcb+a3f46e7fff,gabf8522325+735880ea63,gac2eed3f23+2efcf17695,gb89ab40317+2efcf17695,gbf99507273+8c5ae1fdc5,gd8ff7fe66e+b70e05ef52,gdab6d2f7ff+317eb4c7d4,gdc713202bf+b70e05ef52,gdfd2d52018+b10e285e0f,ge365c994fd+310e8507c4,ge410e46f29+2efcf17695,geaed405ab2+562b3308c0,gffca2db377+8c5ae1fdc5,w.2025.35
LSST Data Management Base Package
Loading...
Searching...
No Matches
repair.py
Go to the documentation of this file.
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__ = ["RepairConfig", "RepairTask", "TooManyCosmicRays"]
23
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
34
35
36class TooManyCosmicRays(pipeBase.AlgorithmError):
37 """Raised if the cosmic ray task fails with too many cosmics.
38
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
50
51 def __str__(self):
52 # Exception doesn't handle **kwargs, so we need a custom str.
53 return f"{self.msg}: {self.metadata}"
54
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
61
62
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 )
82
83 def setDefaults(self):
84 self.interp.useFallbackValueAtEdge = True
85 self.interp.fallbackValueType = "MEANCLIP"
86 self.interp.negativeFallbackAllowed = True
87
88
89class RepairTask(pipeBase.Task):
90 """Repair an exposures defects and cosmic rays via interpolation.
91
92 This task operates on an lsst.afw.image.Exposure in place to
93 interpolate over a set of `~lsst.meas.algorithms.Defect` objects.
94
95 It will also, optionally, find and interpolate any cosmic rays in the lsst.afw.image.Exposure.
96
97 Notes
98 -----
99 Debugging:
100 The available debug variables in RepairTask are:
101
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.
110
111 To investigate the pipe_tasks_repair_Debug, put something like
112
113 .. code-block :: none
114
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
122
123 lsstDebug.Info = DebugInfo
124 into your debug.py file and run runRepair.py with the --debug flag.
125
126 Conversion notes:
127 Display code should be updated once we settle on a standard way of controlling what is displayed.
128 """
129
130 ConfigClass = RepairConfig
131 _DefaultName = "repair"
132
133 def __init__(self, **kwargs):
134 pipeBase.Task.__init__(self, **kwargs)
135 if self.config.doInterpolate:
136 self.makeSubtask("interp")
137
138 @timeMethod
139 def run(self, exposure, defects=None, keepCRs=None):
140 """Repair an Exposure's defects and cosmic rays.
141
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`).
151
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"
164
165 frame = getDebugFrame(self._display, "repair.before")
166 if frame:
167 afwDisplay.Display(frame).mtv(exposure)
168
169 if defects is not None and self.config.doInterpolate:
170 self.interp.run(exposure, defects=defects)
171
172 if self.config.doCosmicRay:
173 self.cosmicRay(exposure, keepCRs=keepCRs)
174
175 frame = getDebugFrame(self._display, "repair.after")
176 if frame:
177 afwDisplay.Display(frame).mtv(exposure)
178
179 def cosmicRay(self, exposure, keepCRs=None):
180 """Mask cosmic rays.
181
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
192
193 assert exposure, "No exposure provided"
194 psf = exposure.getPsf()
195 assert psf, "No psf provided"
196
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
204
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
221
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
235
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())
243
244 except Exception:
245 if display:
246 afwDisplay.Display().mtv(exposure0, title="Failed CR")
247 raise
248
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)
255
256 if display and displayCR:
257 disp = afwDisplay.Display()
258 disp.incrDefaultFrame()
259 disp.mtv(exposure0, title="Post-CR")
260
261 with disp.Buffering():
262 for cr in crs:
263 afwDisplay.utils.drawBBox(cr.getBBox(), borderWidth=0.55)
264
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
cosmicRay(self, exposure, keepCRs=None)
Definition repair.py:179
run(self, exposure, defects=None, keepCRs=None)
Definition repair.py:139
__init__(self, maxCosmicRays, **kwargs)
Definition repair.py:44
Statistics makeStatistics(lsst::afw::image::Image< Pixel > const &img, lsst::afw::image::Mask< image::MaskPixel > const &msk, int const flags, StatisticsControl const &sctrl=StatisticsControl())
Handle a watered-down front-end to the constructor (no variance)
Definition Statistics.h:361