LSST Applications  21.0.0+04719a4bac,21.0.0-1-ga51b5d4+f5e6047307,21.0.0-11-g2b59f77+a9c1acf22d,21.0.0-11-ga42c5b2+86977b0b17,21.0.0-12-gf4ce030+76814010d2,21.0.0-13-g1721dae+760e7a6536,21.0.0-13-g3a573fe+768d78a30a,21.0.0-15-g5a7caf0+f21cbc5713,21.0.0-16-g0fb55c1+b60e2d390c,21.0.0-19-g4cded4ca+71a93a33c0,21.0.0-2-g103fe59+bb20972958,21.0.0-2-g45278ab+04719a4bac,21.0.0-2-g5242d73+3ad5d60fb1,21.0.0-2-g7f82c8f+8babb168e8,21.0.0-2-g8f08a60+06509c8b61,21.0.0-2-g8faa9b5+616205b9df,21.0.0-2-ga326454+8babb168e8,21.0.0-2-gde069b7+5e4aea9c2f,21.0.0-2-gecfae73+1d3a86e577,21.0.0-2-gfc62afb+3ad5d60fb1,21.0.0-25-g1d57be3cd+e73869a214,21.0.0-3-g357aad2+ed88757d29,21.0.0-3-g4a4ce7f+3ad5d60fb1,21.0.0-3-g4be5c26+3ad5d60fb1,21.0.0-3-g65f322c+e0b24896a3,21.0.0-3-g7d9da8d+616205b9df,21.0.0-3-ge02ed75+a9c1acf22d,21.0.0-4-g591bb35+a9c1acf22d,21.0.0-4-g65b4814+b60e2d390c,21.0.0-4-gccdca77+0de219a2bc,21.0.0-4-ge8a399c+6c55c39e83,21.0.0-5-gd00fb1e+05fce91b99,21.0.0-6-gc675373+3ad5d60fb1,21.0.0-64-g1122c245+4fb2b8f86e,21.0.0-7-g04766d7+cd19d05db2,21.0.0-7-gdf92d54+04719a4bac,21.0.0-8-g5674e7b+d1bd76f71f,master-gac4afde19b+a9c1acf22d,w.2021.13
LSST Data Management Base Package
subtractBackground.py
Go to the documentation of this file.
1 #
2 # LSST Data Management System
3 #
4 # Copyright 2008-2017 AURA/LSST.
5 #
6 # This product includes software developed by the
7 # LSST Project (http://www.lsst.org/).
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 LSST License Statement and
20 # the GNU General Public License along with this program. If not,
21 # see <https://www.lsstcorp.org/LegalNotices/>.
22 #
23 __all__ = ("SubtractBackgroundConfig", "SubtractBackgroundTask")
24 
25 import itertools
26 
27 import numpy
28 
29 from lsstDebug import getDebugFrame
30 from lsst.utils import suppress_deprecations
31 
32 import lsst.afw.display as afwDisplay
33 import lsst.afw.image as afwImage
34 import lsst.afw.math as afwMath
35 import lsst.pex.config as pexConfig
36 import lsst.pipe.base as pipeBase
37 
38 
39 class SubtractBackgroundConfig(pexConfig.Config):
40  """!Config for SubtractBackgroundTask
41 
42  @note Many of these fields match fields in lsst.afw.math.BackgroundControl,
43  the control class for lsst.afw.math.makeBackground
44  """
45  statisticsProperty = pexConfig.ChoiceField(
46  doc="type of statistic to use for grid points",
47  dtype=str, default="MEANCLIP",
48  allowed={
49  "MEANCLIP": "clipped mean",
50  "MEAN": "unclipped mean",
51  "MEDIAN": "median",
52  }
53  )
54  undersampleStyle = pexConfig.ChoiceField(
55  doc="behaviour if there are too few points in grid for requested interpolation style",
56  dtype=str, default="REDUCE_INTERP_ORDER",
57  allowed={
58  "THROW_EXCEPTION": "throw an exception if there are too few points",
59  "REDUCE_INTERP_ORDER": "use an interpolation style with a lower order.",
60  "INCREASE_NXNYSAMPLE": "Increase the number of samples used to make the interpolation grid.",
61  },
62  )
63  binSize = pexConfig.RangeField(
64  doc="how large a region of the sky should be used for each background point",
65  dtype=int, default=128, min=1,
66  )
67  binSizeX = pexConfig.RangeField(
68  doc=("Sky region size to be used for each background point in X direction. "
69  "If 0, the binSize config is used."),
70  dtype=int, default=0, min=0,
71  )
72  binSizeY = pexConfig.RangeField(
73  doc=("Sky region size to be used for each background point in Y direction. "
74  "If 0, the binSize config is used."),
75  dtype=int, default=0, min=0,
76  )
77  algorithm = pexConfig.ChoiceField(
78  doc="how to interpolate the background values. This maps to an enum; see afw::math::Background",
79  dtype=str, default="AKIMA_SPLINE", optional=True,
80  allowed={
81  "CONSTANT": "Use a single constant value",
82  "LINEAR": "Use linear interpolation",
83  "NATURAL_SPLINE": "cubic spline with zero second derivative at endpoints",
84  "AKIMA_SPLINE": "higher-level nonlinear spline that is more robust to outliers",
85  "NONE": "No background estimation is to be attempted",
86  },
87  )
88  ignoredPixelMask = pexConfig.ListField(
89  doc="Names of mask planes to ignore while estimating the background",
90  dtype=str, default=["BAD", "EDGE", "DETECTED", "DETECTED_NEGATIVE", "NO_DATA", ],
91  itemCheck=lambda x: x in afwImage.Mask().getMaskPlaneDict().keys(),
92  )
93  isNanSafe = pexConfig.Field(
94  doc="Ignore NaNs when estimating the background",
95  dtype=bool, default=False,
96  )
97 
98  useApprox = pexConfig.Field(
99  doc="Use Approximate (Chebyshev) to model background.",
100  dtype=bool, default=True,
101  )
102  approxOrderX = pexConfig.Field(
103  doc="Approximation order in X for background Chebyshev (valid only with useApprox=True)",
104  dtype=int, default=6,
105  )
106  # Note: Currently X- and Y-orders must be equal due to a limitation in math::Chebyshev1Function2
107  # The following is being added so that the weighting attribute can also be configurable for the
108  # call to afwMath.ApproximateControl
109  approxOrderY = pexConfig.Field(
110  doc="Approximation order in Y for background Chebyshev (valid only with useApprox=True)",
111  dtype=int, default=-1,
112  )
113  weighting = pexConfig.Field(
114  doc="Use inverse variance weighting in calculation (valid only with useApprox=True)",
115  dtype=bool, default=True,
116  )
117 
118 
119 
125 
126 class SubtractBackgroundTask(pipeBase.Task):
127  r"""!Subtract the background from an exposure
128 
129  @anchor SubtractBackgroundTask_
130 
131  @section meas_algorithms_subtractBackground_Contents Contents
132 
133  - @ref meas_algorithms_subtractBackground_Purpose
134  - @ref meas_algorithms_subtractBackground_Initialize
135  - @ref meas_algorithms_subtractBackground_IO
136  - @ref meas_algorithms_subtractBackground_Config
137  - @ref meas_algorithms_subtractBackground_Metadata
138  - @ref meas_algorithms_subtractBackground_Debug
139  - @ref meas_algorithms_subtractBackground_Example
140 
141  @section meas_algorithms_subtractBackground_Purpose Description
142 
143  Fit a model of the background of an exposure and subtract it.
144 
145  @section meas_algorithms_subtractBackground_Initialize Task initialisation
146 
147  @copydoc \_\_init\_\_
148 
149  @section meas_algorithms_subtractBackground_IO Invoking the Task
150 
151  Call `run` to fit the background and subtract it.
152 
153  Call `fitBackground` to fit the background without subtracting it.
154 
155  @section meas_algorithms_subtractBackground_Config Configuration parameters
156 
157  See @ref SubtractBackgroundConfig
158 
159  @section meas_algorithms_subtractBackground_Metadata Quantities set in exposure Metadata
160 
161  The `run` method will optionally set the following items of exposure metadata;
162  the names may be overridden; the defaults are shown:
163  <dl>
164  <dt>BGMEAN <dd>mean value of background
165  <dt>BGVAR <dd>standard deviation of background
166  </dl>
167 
168  @section meas_algorithms_subtractBackground_Debug Debug variables
169 
170  The @link lsst.pipe.base.cmdLineTask.CmdLineTask command line task@endlink interface supports a flag
171  `--debug` to import `debug.py` from your `$PYTHONPATH`; see @ref baseDebug for more about `debug.py`.
172 
173  SubtractBackgroundTask has a debug dictionary containing three integer keys:
174  <dl>
175  <dt>unsubtracted
176  <dd>If >0: `fitBackground` displays the unsubtracted masked image overlaid with the grid of cells
177  used to fit the background in the specified frame
178  <dt>subtracted
179  <dd>If >0: `run` displays the background-subtracted exposure is the specified frame
180  <dt>background
181  <dd>If >0: `run` displays the background image in the specified frame
182  </dl>
183 
184  For example, put something like:
185  @code{.py}
186  import lsstDebug
187  def DebugInfo(name):
188  di = lsstDebug.getInfo(name) # N.b. lsstDebug.Info(name) would call us recursively
189  if name == "lsst.meas.algorithms.subtractBackground":
190  di.display = dict(
191  unsubtracted = 1,
192  subtracted = 2,
193  background = 3,
194  )
195 
196  return di
197 
198  lsstDebug.Info = DebugInfo
199  @endcode
200  into your `debug.py` file and run your task with the `--debug` flag.
201 
202  @section meas_algorithms_subtractBackground_Example A complete example of using SubtractBackgroundTask
203 
204  This code is in @link subtractBackgroundExample.py@endlink in the examples directory, and can be run as:
205  @code
206  python examples/subtractBackgroundExample.py
207  @endcode
208  @dontinclude subtractBackgroundExample.py
209 
210  Import the task (there are some other standard imports; read the file if you're curious)
211  @skipline import SubtractBackgroundTask
212 
213  Create the task, run it, and report mean and variance of background.
214  @skip create the task
215  @until print
216  """
217  ConfigClass = SubtractBackgroundConfig
218  _DefaultName = "subtractBackground"
219 
220  def run(self, exposure, background=None, stats=True, statsKeys=None):
221  """!Fit and subtract the background of an exposure
222 
223  @param[in,out] exposure exposure whose background is to be subtracted
224  @param[in,out] background initial background model already subtracted from exposure
225  (an lsst.afw.math.BackgroundList). May be None if no background has been subtracted.
226  @param[in] stats if True then measure the mean and variance of the full background model
227  and record the results in the exposure's metadata
228  @param[in] statsKeys key names used to store the mean and variance of the background
229  in the exposure's metadata (a pair of strings); if None then use ("BGMEAN", "BGVAR");
230  ignored if stats is false
231 
232  @return an lsst.pipe.base.Struct containing:
233  - background full background model (initial model with changes), an lsst.afw.math.BackgroundList
234  """
235  if background is None:
236  background = afwMath.BackgroundList()
237 
238  maskedImage = exposure.getMaskedImage()
239  fitBg = self.fitBackgroundfitBackground(maskedImage)
240  maskedImage -= fitBg.getImageF(self.config.algorithm, self.config.undersampleStyle)
241 
242  actrl = fitBg.getBackgroundControl().getApproximateControl()
243  background.append((fitBg, getattr(afwMath.Interpolate, self.config.algorithm),
244  fitBg.getAsUsedUndersampleStyle(), actrl.getStyle(),
245  actrl.getOrderX(), actrl.getOrderY(), actrl.getWeighting()))
246 
247  if stats:
248  self._addStats_addStats(exposure, background, statsKeys=statsKeys)
249 
250  subFrame = getDebugFrame(self._display, "subtracted")
251  if subFrame:
252  subDisp = afwDisplay.getDisplay(frame=subFrame)
253  subDisp.mtv(exposure, title="subtracted")
254 
255  bgFrame = getDebugFrame(self._display, "background")
256  if bgFrame:
257  bgDisp = afwDisplay.getDisplay(frame=bgFrame)
258  bgImage = background.getImage()
259  bgDisp.mtv(bgImage, title="background")
260 
261  return pipeBase.Struct(
262  background=background,
263  )
264 
265  def _addStats(self, exposure, background, statsKeys=None):
266  """Add statistics about the background to the exposure's metadata
267 
268  @param[in,out] exposure exposure whose background was subtracted
269  @param[in,out] background background model (an lsst.afw.math.BackgroundList)
270  @param[in] statsKeys key names used to store the mean and variance of the background
271  in the exposure's metadata (a pair of strings); if None then use ("BGMEAN", "BGVAR");
272  ignored if stats is false
273  """
274  netBgImg = background.getImage()
275  if statsKeys is None:
276  statsKeys = ("BGMEAN", "BGVAR")
277  mnkey, varkey = statsKeys
278  meta = exposure.getMetadata()
279  s = afwMath.makeStatistics(netBgImg, afwMath.MEAN | afwMath.VARIANCE)
280  bgmean = s.getValue(afwMath.MEAN)
281  bgvar = s.getValue(afwMath.VARIANCE)
282  meta.addDouble(mnkey, bgmean)
283  meta.addDouble(varkey, bgvar)
284 
285  def fitBackground(self, maskedImage, nx=0, ny=0, algorithm=None):
286  """!Estimate the background of a masked image
287 
288  @param[in] maskedImage masked image whose background is to be computed
289  @param[in] nx number of x bands; if 0 compute from width and config.binSizeX
290  @param[in] ny number of y bands; if 0 compute from height and config.binSizeY
291  @param[in] algorithm name of interpolation algorithm; if None use self.config.algorithm
292 
293  @return fit background as an lsst.afw.math.Background
294 
295  @throw RuntimeError if lsst.afw.math.makeBackground returns None,
296  which is apparently one way it indicates failure
297  """
298 
299  binSizeX = self.config.binSize if self.config.binSizeX == 0 else self.config.binSizeX
300  binSizeY = self.config.binSize if self.config.binSizeY == 0 else self.config.binSizeY
301 
302  if not nx:
303  nx = maskedImage.getWidth()//binSizeX + 1
304  if not ny:
305  ny = maskedImage.getHeight()//binSizeY + 1
306 
307  unsubFrame = getDebugFrame(self._display, "unsubtracted")
308  if unsubFrame:
309  unsubDisp = afwDisplay.getDisplay(frame=unsubFrame)
310  unsubDisp.mtv(maskedImage, title="unsubtracted")
311  xPosts = numpy.rint(numpy.linspace(0, maskedImage.getWidth() + 1, num=nx, endpoint=True))
312  yPosts = numpy.rint(numpy.linspace(0, maskedImage.getHeight() + 1, num=ny, endpoint=True))
313  with unsubDisp.Buffering():
314  for (xMin, xMax), (yMin, yMax) in itertools.product(zip(xPosts[:-1], xPosts[1:]),
315  zip(yPosts[:-1], yPosts[1:])):
316  unsubDisp.line([(xMin, yMin), (xMin, yMax), (xMax, yMax), (xMax, yMin), (xMin, yMin)])
317 
318  sctrl = afwMath.StatisticsControl()
319  badMask = maskedImage.mask.getPlaneBitMask(self.config.ignoredPixelMask)
320 
321  sctrl.setAndMask(badMask)
322  sctrl.setNanSafe(self.config.isNanSafe)
323 
324  self.log.debug("Ignoring mask planes: %s" % ", ".join(self.config.ignoredPixelMask))
325  if (maskedImage.mask.getArray() & badMask).all():
326  raise pipeBase.TaskError("All pixels masked. Cannot estimate background")
327 
328  if algorithm is None:
329  algorithm = self.config.algorithm
330 
331  # TODO: DM-22814. This call to a deprecated BackgroundControl constructor
332  # is necessary to support the algorithm parameter; it # should be replaced with
333  #
334  # afwMath.BackgroundControl(nx, ny, sctrl, self.config.statisticsProperty)
335  #
336  # when algorithm has been deprecated and removed.
337  with suppress_deprecations():
338  bctrl = afwMath.BackgroundControl(algorithm, nx, ny,
339  self.config.undersampleStyle, sctrl,
340  self.config.statisticsProperty)
341 
342  # TODO: The following check should really be done within lsst.afw.math.
343  # With the current code structure, it would need to be accounted for in the doGetImage()
344  # function in BackgroundMI.cc (which currently only checks against the interpolation settings,
345  # which is not appropriate when useApprox=True)
346  # and/or the makeApproximate() function in afw/Approximate.cc.
347  # See ticket DM-2920: "Clean up code in afw for Approximate background
348  # estimation" (which includes a note to remove the following and the
349  # similar checks in pipe_tasks/matchBackgrounds.py once implemented)
350  #
351  # Check that config setting of approxOrder/binSize make sense
352  # (i.e. ngrid (= shortDimension/binSize) > approxOrderX) and perform
353  # appropriate undersampleStlye behavior.
354  if self.config.useApprox:
355  if self.config.approxOrderY not in (self.config.approxOrderX, -1):
356  raise ValueError("Error: approxOrderY not in (approxOrderX, -1)")
357  order = self.config.approxOrderX
358  minNumberGridPoints = order + 1
359  if min(nx, ny) <= order:
360  self.log.warn("Too few points in grid to constrain fit: min(nx, ny) < approxOrder) "
361  "[min(%d, %d) < %d]" % (nx, ny, order))
362  if self.config.undersampleStyle == "THROW_EXCEPTION":
363  raise ValueError("Too few points in grid (%d, %d) for order (%d) and binSize (%d, %d)" %
364  (nx, ny, order, binSizeX, binSizeY))
365  elif self.config.undersampleStyle == "REDUCE_INTERP_ORDER":
366  if order < 1:
367  raise ValueError("Cannot reduce approxOrder below 0. "
368  "Try using undersampleStyle = \"INCREASE_NXNYSAMPLE\" instead?")
369  order = min(nx, ny) - 1
370  self.log.warn("Reducing approxOrder to %d" % order)
371  elif self.config.undersampleStyle == "INCREASE_NXNYSAMPLE":
372  # Reduce bin size to the largest acceptable square bins
373  newBinSize = min(maskedImage.getWidth(), maskedImage.getHeight())//(minNumberGridPoints-1)
374  if newBinSize < 1:
375  raise ValueError("Binsize must be greater than 0")
376  newNx = maskedImage.getWidth()//newBinSize + 1
377  newNy = maskedImage.getHeight()//newBinSize + 1
378  bctrl.setNxSample(newNx)
379  bctrl.setNySample(newNy)
380  self.log.warn("Decreasing binSize from (%d, %d) to %d for a grid of (%d, %d)" %
381  (binSizeX, binSizeY, newBinSize, newNx, newNy))
382 
383  actrl = afwMath.ApproximateControl(afwMath.ApproximateControl.CHEBYSHEV, order, order,
384  self.config.weighting)
385  bctrl.setApproximateControl(actrl)
386 
387  bg = afwMath.makeBackground(maskedImage, bctrl)
388  if bg is None:
389  raise RuntimeError("lsst.afw.math.makeBackground failed to fit a background model")
390  return bg
int min
Represent a 2-dimensional array of bitmask pixels.
Definition: Mask.h:77
Control how to make an approximation.
Definition: Approximate.h:48
Pass parameters to a Background object.
Definition: Background.h:56
Pass parameters to a Statistics object.
Definition: Statistics.h:93
def fitBackground(self, maskedImage, nx=0, ny=0, algorithm=None)
Estimate the background of a masked image.
def run(self, exposure, background=None, stats=True, statsKeys=None)
Fit and subtract the background of an exposure.
def _addStats(self, exposure, background, statsKeys=None)
Backwards-compatibility support for depersisting the old Calib (FluxMag0/FluxMag0Err) objects.
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:354
std::shared_ptr< Background > makeBackground(ImageT const &img, BackgroundControl const &bgCtrl)
A convenience function that uses function overloading to make the correct type of Background.
Definition: Background.h:526
bool all(CoordinateExpr< N > const &expr) noexcept
Return true if all elements are true.
def suppress_deprecations(category=FutureWarning)
Definition: deprecated.py:79
def getDebugFrame(debugDisplay, name)
Definition: lsstDebug.py:95