lsst.ip.diffim  16.0-8-g4aca173+17
dcrModel.py
Go to the documentation of this file.
1 # This file is part of ip_diffim.
2 #
3 # LSST Data Management System
4 # This product includes software developed by the
5 # LSST Project (http://www.lsst.org/).
6 # See COPYRIGHT file at the top of the source tree.
7 #
8 # This program is free software: you can redistribute it and/or modify
9 # it under the terms of the GNU General Public License as published by
10 # the Free Software Foundation, either version 3 of the License, or
11 # (at your option) any later version.
12 #
13 # This program is distributed in the hope that it will be useful,
14 # but WITHOUT ANY WARRANTY; without even the implied warranty of
15 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 # GNU General Public License for more details.
17 #
18 # You should have received a copy of the LSST License Statement and
19 # the GNU General Public License along with this program. If not,
20 # see <https://www.lsstcorp.org/LegalNotices/>.
21 #
22 
23 import numpy as np
24 from lsst.afw.coord.refraction import differentialRefraction
25 import lsst.afw.geom as afwGeom
26 from lsst.afw.geom import AffineTransform
27 from lsst.afw.geom import makeTransform
28 import lsst.afw.image as afwImage
29 import lsst.afw.math as afwMath
30 from lsst.geom import radians
31 
32 __all__ = ["DcrModel", "applyDcr", "calculateDcr", "calculateImageParallacticAngle"]
33 
34 
35 class DcrModel:
36  """A model of the true sky after correcting chromatic effects.
37 
38  Attributes
39  ----------
40  dcrNumSubfilters : `int`
41  Number of sub-filters used to model chromatic effects within a band.
42  filterInfo : `lsst.afw.image.Filter`
43  The filter definition, set in the current instruments' obs package.
44  modelImages : `list` of `lsst.afw.image.MaskedImage`
45  A list of masked images, each containing the model for one subfilter
46 
47  Parameters
48  ----------
49  modelImages : `list` of `lsst.afw.image.MaskedImage`
50  A list of masked images, each containing the model for one subfilter.
51  filterInfo : `lsst.afw.image.Filter`, optional
52  The filter definition, set in the current instruments' obs package.
53  Required for any calculation of DCR, including making matched templates.
54 
55  Notes
56  -----
57  The ``DcrModel`` contains an estimate of the true sky, at a higher
58  wavelength resolution than the input observations. It can be forward-
59  modeled to produce Differential Chromatic Refraction (DCR) matched
60  templates for a given ``Exposure``, and provides utilities for conditioning
61  the model in ``dcrAssembleCoadd`` to avoid oscillating solutions between
62  iterations of forward modeling or between the subfilters of the model.
63  """
64 
65  def __init__(self, modelImages, filterInfo=None, psf=None):
66  self.dcrNumSubfilters = len(modelImages)
67  self.modelImages = modelImages
68  self._filter = filterInfo
69  self._psf = psf
70 
71  @classmethod
72  def fromImage(cls, maskedImage, dcrNumSubfilters, filterInfo=None, psf=None):
73  """Initialize a DcrModel by dividing a coadd between the subfilters.
74 
75  Parameters
76  ----------
77  maskedImage : `lsst.afw.image.MaskedImage`
78  Input coadded image to divide equally between the subfilters.
79  dcrNumSubfilters : `int`
80  Number of sub-filters used to model chromatic effects within a band.
81  filterInfo : `lsst.afw.image.Filter`, optional
82  The filter definition, set in the current instruments' obs package.
83  Required for any calculation of DCR, including making matched templates.
84  psf : `lsst.afw.detection.Psf`, optional
85  Point spread function (PSF) of the model.
86  Required if the ``DcrModel`` will be persisted.
87 
88  Returns
89  -------
90  dcrModel : `lsst.pipe.tasks.DcrModel`
91  Best fit model of the true sky after correcting chromatic effects.
92  """
93  # NANs will potentially contaminate the entire image,
94  # depending on the shift or convolution type used.
95  model = maskedImage.clone()
96  badPixels = np.isnan(model.image.array) | np.isnan(model.variance.array)
97  model.image.array[badPixels] = 0.
98  model.variance.array[badPixels] = 0.
99  model.image.array /= dcrNumSubfilters
100  # We divide the variance by N and not N**2 because we will assume each
101  # subfilter is independent. That means that the significance of
102  # detected sources will be lower by a factor of sqrt(N) in the
103  # subfilter images, but we will recover it when we combine the
104  # subfilter images to construct matched templates.
105  model.variance.array /= dcrNumSubfilters
106  model.mask.array[badPixels] = model.mask.getPlaneBitMask("NO_DATA")
107  modelImages = [model, ]
108  for subfilter in range(1, dcrNumSubfilters):
109  modelImages.append(model.clone())
110  return cls(modelImages, filterInfo, psf)
111 
112  @classmethod
113  def fromDataRef(cls, dataRef, datasetType="dcrCoadd", numSubfilters=None, **kwargs):
114  """Load an existing DcrModel from a repository.
115 
116  Parameters
117  ----------
118  dataRef : `lsst.daf.persistence.ButlerDataRef`
119  Data reference defining the patch for coaddition and the
120  reference Warp
121  datasetType : `str`, optional
122  Name of the DcrModel in the registry {"dcrCoadd", "dcrCoadd_sub"}
123  numSubfilters : `int`
124  Number of sub-filters used to model chromatic effects within a band.
125  **kwargs
126  Additional keyword arguments to pass to look up the model in the data registry.
127  Common keywords and their types include: ``tract``:`str`, ``patch``:`str`,
128  ``bbox``:`lsst.afw.geom.Box2I`
129 
130  Returns
131  -------
132  dcrModel : `lsst.pipe.tasks.DcrModel`
133  Best fit model of the true sky after correcting chromatic effects.
134  """
135  modelImages = []
136  filterInfo = None
137  psf = None
138  for subfilter in range(numSubfilters):
139  dcrCoadd = dataRef.get(datasetType, subfilter=subfilter,
140  numSubfilters=numSubfilters, **kwargs)
141  if filterInfo is None:
142  filterInfo = dcrCoadd.getFilter()
143  if psf is None:
144  psf = dcrCoadd.getPsf()
145  modelImages.append(dcrCoadd.maskedImage)
146  return cls(modelImages, filterInfo, psf)
147 
148  def __len__(self):
149  """Return the number of subfilters.
150 
151  Returns
152  -------
153  dcrNumSubfilters : `int`
154  The number of DCR subfilters in the model.
155  """
156  return self.dcrNumSubfilters
157 
158  def __getitem__(self, subfilter):
159  """Iterate over the subfilters of the DCR model.
160 
161  Parameters
162  ----------
163  subfilter : `int`
164  Index of the current ``subfilter`` within the full band.
165  Negative indices are allowed, and count in reverse order
166  from the highest ``subfilter``.
167 
168  Returns
169  -------
170  modelImage : `lsst.afw.image.MaskedImage`
171  The DCR model for the given ``subfilter``.
172 
173  Raises
174  ------
175  IndexError
176  If the requested ``subfilter`` is greater or equal to the number
177  of subfilters in the model.
178  """
179  if np.abs(subfilter) >= len(self):
180  raise IndexError("subfilter out of bounds.")
181  return self.modelImages[subfilter]
182 
183  def __setitem__(self, subfilter, maskedImage):
184  """Update the model image for one subfilter.
185 
186  Parameters
187  ----------
188  subfilter : `int`
189  Index of the current subfilter within the full band.
190  maskedImage : `lsst.afw.image.MaskedImage`
191  The DCR model to set for the given ``subfilter``.
192 
193  Raises
194  ------
195  IndexError
196  If the requested ``subfilter`` is greater or equal to the number
197  of subfilters in the model.
198  ValueError
199  If the bounding box of the new image does not match.
200  """
201  if np.abs(subfilter) >= len(self):
202  raise IndexError("subfilter out of bounds.")
203  if maskedImage.getBBox() != self.getBBox():
204  raise ValueError("The bounding box of a subfilter must not change.")
205  self.modelImages[subfilter] = maskedImage
206 
207  @property
208  def filter(self):
209  """Return the filter of the model.
210 
211  Returns
212  -------
213  filter : `lsst.afw.image.Filter`
214  The filter definition, set in the current instruments' obs package.
215  """
216  return self._filter
217 
218  @property
219  def psf(self):
220  """Return the psf of the model.
221 
222  Returns
223  -------
224  psf : `lsst.afw.detection.Psf`
225  Point spread function (PSF) of the model.
226  """
227  return self._psf
228 
229  @property
230  def bbox(self):
231  """Return the common bounding box of each subfilter image.
232 
233  Returns
234  -------
235  bbox : `lsst.afw.geom.Box2I`
236  Bounding box of the DCR model.
237  """
238  return self[0].getBBox()
239 
240  def getReferenceImage(self, bbox=None):
241  """Create a simple template from the DCR model.
242 
243  Parameters
244  ----------
245  bbox : `lsst.afw.geom.Box2I`, optional
246  Sub-region of the coadd. Returns the entire image if `None`.
247 
248  Returns
249  -------
250  templateImage : `numpy.ndarray`
251  The template with no chromatic effects applied.
252  """
253  return np.mean([model[bbox].image.array for model in self], axis=0)
254 
255  def assign(self, dcrSubModel, bbox=None):
256  """Update a sub-region of the ``DcrModel`` with new values.
257 
258  Parameters
259  ----------
260  dcrSubModel : `lsst.pipe.tasks.DcrModel`
261  New model of the true scene after correcting chromatic effects.
262  bbox : `lsst.afw.geom.Box2I`, optional
263  Sub-region of the coadd.
264  Defaults to the bounding box of ``dcrSubModel``.
265 
266  Raises
267  ------
268  ValueError
269  If the new model has a different number of subfilters.
270  """
271  if len(dcrSubModel) != len(self):
272  raise ValueError("The number of DCR subfilters must be the same "
273  "between the old and new models.")
274  if bbox is None:
275  bbox = self.getBBox()
276  for model, subModel in zip(self, dcrSubModel):
277  model.assign(subModel[bbox], bbox)
278 
279  def buildMatchedTemplate(self, exposure=None, warpCtrl=None,
280  visitInfo=None, bbox=None, wcs=None, mask=None):
281  """Create a DCR-matched template image for an exposure.
282 
283  Parameters
284  ----------
285  exposure : `lsst.afw.image.Exposure`, optional
286  The input exposure to build a matched template for.
287  May be omitted if all of the metadata is supplied separately
288  warpCtrl : `lsst.afw.Math.WarpingControl`, optional
289  Configuration settings for warping an image.
290  If not set, defaults to a lanczos3 warping kernel for the image,
291  and a bilinear kernel for the mask
292  visitInfo : `lsst.afw.image.VisitInfo`, optional
293  Metadata for the exposure. Ignored if ``exposure`` is set.
294  bbox : `lsst.afw.geom.Box2I`, optional
295  Sub-region of the coadd. Ignored if ``exposure`` is set.
296  wcs : `lsst.afw.geom.SkyWcs`, optional
297  Coordinate system definition (wcs) for the exposure.
298  Ignored if ``exposure`` is set.
299  mask : `lsst.afw.image.Mask`, optional
300  reference mask to use for the template image.
301 
302  Returns
303  -------
304  templateImage : `lsst.afw.image.maskedImageF`
305  The DCR-matched template
306 
307  Raises
308  ------
309  ValueError
310  If neither ``exposure`` or all of ``visitInfo``, ``bbox``, and ``wcs`` are set.
311  """
312  if self.filter is None:
313  raise ValueError("'filterInfo' must be set for the DcrModel in order to calculate DCR.")
314  if exposure is not None:
315  visitInfo = exposure.getInfo().getVisitInfo()
316  bbox = exposure.getBBox()
317  wcs = exposure.getInfo().getWcs()
318  elif visitInfo is None or bbox is None or wcs is None:
319  raise ValueError("Either exposure or visitInfo, bbox, and wcs must be set.")
320  if warpCtrl is None:
321  # Turn off the warping cache, since we set the linear interpolation length to the entire subregion
322  # This warper is only used for applying DCR shifts, which are assumed to be uniform across a patch
323  warpCtrl = afwMath.WarpingControl("lanczos3", "bilinear",
324  cacheSize=0, interpLength=max(bbox.getDimensions()))
325 
326  dcrShift = calculateDcr(visitInfo, wcs, self.filter, len(self))
327  templateImage = afwImage.MaskedImageF(bbox)
328  for subfilter, dcr in enumerate(dcrShift):
329  templateImage += applyDcr(self[subfilter][bbox], dcr, warpCtrl)
330  if mask is not None:
331  templateImage.setMask(mask[bbox])
332  return templateImage
333 
334  def buildMatchedExposure(self, exposure=None, warpCtrl=None,
335  visitInfo=None, bbox=None, wcs=None, mask=None):
336  """Wrapper to create an exposure from a template image.
337 
338  Parameters
339  ----------
340  exposure : `lsst.afw.image.Exposure`, optional
341  The input exposure to build a matched template for.
342  May be omitted if all of the metadata is supplied separately
343  warpCtrl : `lsst.afw.Math.WarpingControl`
344  Configuration settings for warping an image
345  visitInfo : `lsst.afw.image.VisitInfo`, optional
346  Metadata for the exposure. Ignored if ``exposure`` is set.
347  bbox : `lsst.afw.geom.Box2I`, optional
348  Sub-region of the coadd. Ignored if ``exposure`` is set.
349  wcs : `lsst.afw.geom.SkyWcs`, optional
350  Coordinate system definition (wcs) for the exposure.
351  Ignored if ``exposure`` is set.
352  mask : `lsst.afw.image.Mask`, optional
353  reference mask to use for the template image.
354 
355  Returns
356  -------
357  templateExposure : `lsst.afw.image.exposureF`
358  The DCR-matched template
359  """
360  templateImage = self.buildMatchedTemplate(exposure, warpCtrl, visitInfo, bbox, wcs, mask)
361  templateExposure = afwImage.ExposureF(bbox, wcs)
362  templateExposure.setMaskedImage(templateImage)
363  templateExposure.setPsf(self.psf)
364  templateExposure.setFilter(self.filter)
365  return templateExposure
366 
367  def conditionDcrModel(self, subfilter, newModel, bbox, gain=1.):
368  """Average two iterations' solutions to reduce oscillations.
369 
370  Parameters
371  ----------
372  subfilter : `int`
373  Index of the current subfilter within the full band.
374  newModel : `lsst.afw.image.MaskedImage`
375  The new DCR model for one subfilter from the current iteration.
376  Values in ``newModel`` that are extreme compared with the last
377  iteration are modified in place.
378  bbox : `lsst.afw.geom.Box2I`
379  Sub-region of the coadd
380  gain : `float`, optional
381  Additional weight to apply to the model from the current iteration.
382  Defaults to 1.0, which gives equal weight to both solutions.
383  """
384  # Calculate weighted averages of the image and variance planes.
385  # Note that ``newModel *= gain`` would multiply the variance by ``gain**2``
386  newModel.image *= gain
387  newModel.image += self[subfilter][bbox].image
388  newModel.image /= 1. + gain
389  newModel.variance *= gain
390  newModel.variance += self[subfilter][bbox].variance
391  newModel.variance /= 1. + gain
392 
393  def clampModel(self, subfilter, newModel, bbox, statsCtrl, regularizeSigma, modelClampFactor,
394  convergenceMaskPlanes="DETECTED"):
395  """Restrict large variations in the model between iterations.
396 
397  Parameters
398  ----------
399  subfilter : `int`
400  Index of the current subfilter within the full band.
401  newModel : `lsst.afw.image.MaskedImage`
402  The new DCR model for one subfilter from the current iteration.
403  Values in ``newModel`` that are extreme compared with the last
404  iteration are modified in place.
405  bbox : `lsst.afw.geom.Box2I`
406  Sub-region to coadd
407  statsCtrl : `lsst.afw.math.StatisticsControl`
408  Statistics control object for coadd
409  regularizeSigma : `float`
410  Threshold to exclude noise-like pixels from regularization.
411  modelClampFactor : `float`
412  Maximum relative change of the model allowed between iterations.
413  convergenceMaskPlanes : `list` of `str`, or `str`, optional
414  Mask planes to use to calculate convergence.
415  Default value is set in ``calculateNoiseCutoff`` if not supplied.
416  """
417  newImage = newModel.image.array
418  oldImage = self[subfilter][bbox].image.array
419  noiseCutoff = self.calculateNoiseCutoff(newModel, statsCtrl, regularizeSigma,
420  convergenceMaskPlanes=convergenceMaskPlanes)
421  # Catch any invalid values
422  nanPixels = np.isnan(newImage)
423  newImage[nanPixels] = 0.
424  infPixels = np.isinf(newImage)
425  newImage[infPixels] = oldImage[infPixels]*modelClampFactor
426  # Clip pixels that have very high amplitude, compared with the previous iteration.
427  clampPixels = np.abs(newImage - oldImage) > (np.abs(oldImage*(modelClampFactor - 1)) +
428  noiseCutoff)
429  # Set high amplitude pixels to a multiple or fraction of the old model value,
430  # depending on whether the new model is higher or lower than the old
431  highPixels = newImage > oldImage
432  newImage[clampPixels & highPixels] = oldImage[clampPixels & highPixels]*modelClampFactor
433  newImage[clampPixels & ~highPixels] = oldImage[clampPixels & ~highPixels]/modelClampFactor
434 
435  def regularizeModel(self, bbox, mask, statsCtrl, regularizeSigma, clampFrequency,
436  convergenceMaskPlanes="DETECTED"):
437  """Restrict large variations in the model between subfilters.
438 
439  Any flux subtracted by the restriction is accumulated from all
440  subfilters, and divided evenly to each afterwards in order to preserve
441  total flux.
442 
443  Parameters
444  ----------
445  bbox : `lsst.afw.geom.Box2I`
446  Sub-region to coadd
447  mask : `lsst.afw.image.Mask`
448  Reference mask to use for all model planes.
449  statsCtrl : `lsst.afw.math.StatisticsControl`
450  Statistics control object for coadd
451  regularizeSigma : `float`
452  Threshold to exclude noise-like pixels from regularization.
453  clampFrequency : `float`
454  Maximum relative change of the model allowed between subfilters.
455  convergenceMaskPlanes : `list` of `str`, or `str`, optional
456  Mask planes to use to calculate convergence. (Default is "DETECTED")
457  Default value is set in ``calculateNoiseCutoff`` if not supplied.
458  """
459  templateImage = self.getReferenceImage(bbox)
460  excess = np.zeros_like(templateImage)
461  for model in self:
462  noiseCutoff = self.calculateNoiseCutoff(model, statsCtrl, regularizeSigma,
463  convergenceMaskPlanes=convergenceMaskPlanes,
464  mask=mask[bbox])
465  modelVals = model.image.array
466  highPixels = (modelVals > (templateImage*clampFrequency + noiseCutoff))
467  excess[highPixels] += modelVals[highPixels] - templateImage[highPixels]*clampFrequency
468  modelVals[highPixels] = templateImage[highPixels]*clampFrequency
469  lowPixels = (modelVals < templateImage/clampFrequency - noiseCutoff)
470  excess[lowPixels] += modelVals[lowPixels] - templateImage[lowPixels]/clampFrequency
471  modelVals[lowPixels] = templateImage[lowPixels]/clampFrequency
472  excess /= len(self)
473  for model in self:
474  model.image.array += excess
475 
476  def calculateNoiseCutoff(self, maskedImage, statsCtrl, regularizeSigma,
477  convergenceMaskPlanes="DETECTED", mask=None):
478  """Helper function to calculate the background noise level of an image.
479 
480  Parameters
481  ----------
482  maskedImage : `lsst.afw.image.MaskedImage`
483  The input image to evaluate the background noise properties.
484  statsCtrl : `lsst.afw.math.StatisticsControl`
485  Statistics control object for coadd
486  regularizeSigma : `float`
487  Threshold to exclude noise-like pixels from regularization.
488  convergenceMaskPlanes : `list` of `str`, or `str`
489  Mask planes to use to calculate convergence.
490  mask : `lsst.afw.image.Mask`, Optional
491  Optional alternate mask
492 
493  Returns
494  -------
495  noiseCutoff : `float`
496  The threshold value to treat pixels as noise in an image..
497  """
498  convergeMask = maskedImage.mask.getPlaneBitMask(convergenceMaskPlanes)
499  if mask is None:
500  mask = maskedImage.mask
501  backgroundPixels = mask.array & (statsCtrl.getAndMask() | convergeMask) == 0
502  noiseCutoff = regularizeSigma*np.std(maskedImage.image.array[backgroundPixels])
503  return noiseCutoff
504 
505 
506 def applyDcr(maskedImage, dcr, warpCtrl, bbox=None, useInverse=False):
507  """Shift a masked image.
508 
509  Parameters
510  ----------
511  maskedImage : `lsst.afw.image.MaskedImage`
512  The input masked image to shift.
513  dcr : `lsst.afw.geom.Extent2I`
514  Shift calculated with ``calculateDcr``.
515  warpCtrl : `lsst.afw.math.WarpingControl`
516  Configuration settings for warping an image
517  bbox : `lsst.afw.geom.Box2I`, optional
518  Sub-region of the masked image to shift.
519  Shifts the entire image if None (Default).
520  useInverse : `bool`, optional
521  Use the reverse of ``dcr`` for the shift. Default: False
522 
523  Returns
524  -------
525  `lsst.afw.image.maskedImageF`
526  A masked image, with the pixels within the bounding box shifted.
527  """
528  padValue = afwImage.pixel.SinglePixelF(0., maskedImage.mask.getPlaneBitMask("NO_DATA"), 0)
529  if bbox is None:
530  bbox = maskedImage.getBBox()
531  shiftedImage = afwImage.MaskedImageF(bbox)
532  transform = makeTransform(AffineTransform((-1.0 if useInverse else 1.0)*dcr))
533  afwMath.warpImage(shiftedImage, maskedImage[bbox],
534  transform, warpCtrl, padValue=padValue)
535  return shiftedImage
536 
537 
538 def calculateDcr(visitInfo, wcs, filterInfo, dcrNumSubfilters):
539  """Calculate the shift in pixels of an exposure due to DCR.
540 
541  Parameters
542  ----------
543  visitInfo : `lsst.afw.image.VisitInfo`
544  Metadata for the exposure.
545  wcs : `lsst.afw.geom.SkyWcs`
546  Coordinate system definition (wcs) for the exposure.
547  filterInfo : `lsst.afw.image.Filter`
548  The filter definition, set in the current instruments' obs package.
549  dcrNumSubfilters : `int`
550  Number of sub-filters used to model chromatic effects within a band.
551 
552  Returns
553  -------
554  `lsst.afw.geom.Extent2I`
555  The 2D shift due to DCR, in pixels.
556  """
557  rotation = calculateImageParallacticAngle(visitInfo, wcs)
558  dcrShift = []
559  lambdaEff = filterInfo.getFilterProperty().getLambdaEff()
560  for wl0, wl1 in wavelengthGenerator(filterInfo, dcrNumSubfilters):
561  # Note that diffRefractAmp can be negative, since it's relative to the midpoint of the full band
562  diffRefractAmp0 = differentialRefraction(wl0, lambdaEff,
563  elevation=visitInfo.getBoresightAzAlt().getLatitude(),
564  observatory=visitInfo.getObservatory(),
565  weather=visitInfo.getWeather())
566  diffRefractAmp1 = differentialRefraction(wl1, lambdaEff,
567  elevation=visitInfo.getBoresightAzAlt().getLatitude(),
568  observatory=visitInfo.getObservatory(),
569  weather=visitInfo.getWeather())
570  diffRefractAmp = (diffRefractAmp0 + diffRefractAmp1)/2.
571  diffRefractPix = diffRefractAmp.asArcseconds()/wcs.getPixelScale().asArcseconds()
572  dcrShift.append(afwGeom.Extent2D(diffRefractPix*np.cos(rotation.asRadians()),
573  diffRefractPix*np.sin(rotation.asRadians())))
574  return dcrShift
575 
576 
577 def calculateImageParallacticAngle(visitInfo, wcs):
578  """Calculate the total sky rotation angle of an exposure.
579 
580  Parameters
581  ----------
582  visitInfo : `lsst.afw.image.VisitInfo`
583  Metadata for the exposure.
584  wcs : `lsst.afw.geom.SkyWcs`
585  Coordinate system definition (wcs) for the exposure.
586 
587  Returns
588  -------
589  `lsst.geom.Angle`
590  The rotation of the image axis, East from North.
591  Equal to the parallactic angle plus any additional rotation of the
592  coordinate system.
593  A rotation angle of 0 degrees is defined with
594  North along the +y axis and East along the +x axis.
595  A rotation angle of 90 degrees is defined with
596  North along the +x axis and East along the -y axis.
597  """
598  parAngle = visitInfo.getBoresightParAngle().asRadians()
599  cd = wcs.getCdMatrix()
600  cdAngle = (np.arctan2(-cd[0, 1], cd[0, 0]) + np.arctan2(cd[1, 0], cd[1, 1]))/2.
601  rotAngle = (cdAngle + parAngle)*radians
602  return rotAngle
603 
604 
605 def wavelengthGenerator(filterInfo, dcrNumSubfilters):
606  """Iterate over the wavelength endpoints of subfilters.
607 
608  Parameters
609  ----------
610  filterInfo : `lsst.afw.image.Filter`
611  The filter definition, set in the current instruments' obs package.
612  dcrNumSubfilters : `int`
613  Number of sub-filters used to model chromatic effects within a band.
614 
615  Yields
616  ------
617  `tuple` of two `float`
618  The next set of wavelength endpoints for a subfilter, in nm.
619  """
620  lambdaMin = filterInfo.getFilterProperty().getLambdaMin()
621  lambdaMax = filterInfo.getFilterProperty().getLambdaMax()
622  wlStep = (lambdaMax - lambdaMin)/dcrNumSubfilters
623  for wl in np.linspace(lambdaMin, lambdaMax, dcrNumSubfilters, endpoint=False):
624  yield (wl, wl + wlStep)
def fromImage(cls, maskedImage, dcrNumSubfilters, filterInfo=None, psf=None)
Definition: dcrModel.py:72
def __setitem__(self, subfilter, maskedImage)
Definition: dcrModel.py:183
def buildMatchedTemplate(self, exposure=None, warpCtrl=None, visitInfo=None, bbox=None, wcs=None, mask=None)
Definition: dcrModel.py:280
def calculateImageParallacticAngle(visitInfo, wcs)
Definition: dcrModel.py:577
int warpImage(DestImageT &destImage, SrcImageT const &srcImage, geom::TransformPoint2ToPoint2 const &srcToDest, WarpingControl const &control, typename DestImageT::SinglePixel padValue=lsst::afw::math::edgePixel< DestImageT >(typename lsst::afw::image::detail::image_traits< DestImageT >::image_category()))
def calculateNoiseCutoff(self, maskedImage, statsCtrl, regularizeSigma, convergenceMaskPlanes="DETECTED", mask=None)
Definition: dcrModel.py:477
def __getitem__(self, subfilter)
Definition: dcrModel.py:158
def __init__(self, modelImages, filterInfo=None, psf=None)
Definition: dcrModel.py:65
def wavelengthGenerator(filterInfo, dcrNumSubfilters)
Definition: dcrModel.py:605
def buildMatchedExposure(self, exposure=None, warpCtrl=None, visitInfo=None, bbox=None, wcs=None, mask=None)
Definition: dcrModel.py:335
def calculateDcr(visitInfo, wcs, filterInfo, dcrNumSubfilters)
Definition: dcrModel.py:538
def fromDataRef(cls, dataRef, datasetType="dcrCoadd", numSubfilters=None, kwargs)
Definition: dcrModel.py:113
def conditionDcrModel(self, subfilter, newModel, bbox, gain=1.)
Definition: dcrModel.py:367
def getReferenceImage(self, bbox=None)
Definition: dcrModel.py:240
def assign(self, dcrSubModel, bbox=None)
Definition: dcrModel.py:255
def applyDcr(maskedImage, dcr, warpCtrl, bbox=None, useInverse=False)
Definition: dcrModel.py:506
def regularizeModel(self, bbox, mask, statsCtrl, regularizeSigma, clampFrequency, convergenceMaskPlanes="DETECTED")
Definition: dcrModel.py:436
def clampModel(self, subfilter, newModel, bbox, statsCtrl, regularizeSigma, modelClampFactor, convergenceMaskPlanes="DETECTED")
Definition: dcrModel.py:394