LSST Applications  21.0.0-172-gfb10e10a+18fedfabac,22.0.0+297cba6710,22.0.0+80564b0ff1,22.0.0+8d77f4f51a,22.0.0+a28f4c53b1,22.0.0+dcf3732eb2,22.0.1-1-g7d6de66+2a20fdde0d,22.0.1-1-g8e32f31+297cba6710,22.0.1-1-geca5380+7fa3b7d9b6,22.0.1-12-g44dc1dc+2a20fdde0d,22.0.1-15-g6a90155+515f58c32b,22.0.1-16-g9282f48+790f5f2caa,22.0.1-2-g92698f7+dcf3732eb2,22.0.1-2-ga9b0f51+7fa3b7d9b6,22.0.1-2-gd1925c9+bf4f0e694f,22.0.1-24-g1ad7a390+a9625a72a8,22.0.1-25-g5bf6245+3ad8ecd50b,22.0.1-25-gb120d7b+8b5510f75f,22.0.1-27-g97737f7+2a20fdde0d,22.0.1-32-gf62ce7b1+aa4237961e,22.0.1-4-g0b3f228+2a20fdde0d,22.0.1-4-g243d05b+871c1b8305,22.0.1-4-g3a563be+32dcf1063f,22.0.1-4-g44f2e3d+9e4ab0f4fa,22.0.1-42-gca6935d93+ba5e5ca3eb,22.0.1-5-g15c806e+85460ae5f3,22.0.1-5-g58711c4+611d128589,22.0.1-5-g75bb458+99c117b92f,22.0.1-6-g1c63a23+7fa3b7d9b6,22.0.1-6-g50866e6+84ff5a128b,22.0.1-6-g8d3140d+720564cf76,22.0.1-6-gd805d02+cc5644f571,22.0.1-8-ge5750ce+85460ae5f3,master-g6e05de7fdc+babf819c66,master-g99da0e417a+8d77f4f51a,w.2021.48
LSST Data Management Base Package
isrMock.py
Go to the documentation of this file.
1 # This file is part of ip_isr.
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 import copy
23 import numpy as np
24 import tempfile
25 
26 import lsst.geom
27 import lsst.afw.geom as afwGeom
28 import lsst.afw.image as afwImage
29 
30 import lsst.afw.cameraGeom.utils as afwUtils
31 import lsst.afw.cameraGeom.testUtils as afwTestUtils
32 import lsst.pex.config as pexConfig
33 import lsst.pipe.base as pipeBase
34 from .crosstalk import CrosstalkCalib
35 from .defects import Defects
36 
37 __all__ = ["IsrMockConfig", "IsrMock", "RawMock", "TrimmedRawMock", "RawDictMock",
38  "CalibratedRawMock", "MasterMock",
39  "BiasMock", "DarkMock", "FlatMock", "FringeMock", "UntrimmedFringeMock",
40  "BfKernelMock", "DefectMock", "CrosstalkCoeffMock", "TransmissionMock",
41  "DataRefMock"]
42 
43 
44 class IsrMockConfig(pexConfig.Config):
45  """Configuration parameters for isrMock.
46 
47  These parameters produce generic fixed position signals from
48  various sources, and combine them in a way that matches how those
49  signals are combined to create real data. The camera used is the
50  test camera defined by the afwUtils code.
51  """
52  # Detector parameters. "Exposure" parameters.
53  isLsstLike = pexConfig.Field(
54  dtype=bool,
55  default=False,
56  doc="If True, products have one raw image per amplifier, otherwise, one raw image per detector.",
57  )
58  plateScale = pexConfig.Field(
59  dtype=float,
60  default=20.0,
61  doc="Plate scale used in constructing mock camera.",
62  )
63  radialDistortion = pexConfig.Field(
64  dtype=float,
65  default=0.925,
66  doc="Radial distortion term used in constructing mock camera.",
67  )
68  isTrimmed = pexConfig.Field(
69  dtype=bool,
70  default=True,
71  doc="If True, amplifiers have been trimmed and mosaicked to remove regions outside the data BBox.",
72  )
73  detectorIndex = pexConfig.Field(
74  dtype=int,
75  default=20,
76  doc="Index for the detector to use. The default value uses a standard 2x4 array of amps.",
77  )
78  rngSeed = pexConfig.Field(
79  dtype=int,
80  default=20000913,
81  doc="Seed for random number generator used to add noise.",
82  )
83  # TODO: DM-18345 Check that mocks scale correctly when gain != 1.0
84  gain = pexConfig.Field(
85  dtype=float,
86  default=1.0,
87  doc="Gain for simulated data in e^-/DN.",
88  )
89  readNoise = pexConfig.Field(
90  dtype=float,
91  default=5.0,
92  doc="Read noise of the detector in e-.",
93  )
94  expTime = pexConfig.Field(
95  dtype=float,
96  default=5.0,
97  doc="Exposure time for simulated data.",
98  )
99 
100  # Signal parameters
101  skyLevel = pexConfig.Field(
102  dtype=float,
103  default=1000.0,
104  doc="Background contribution to be generated from 'the sky' in DN.",
105  )
106  sourceFlux = pexConfig.ListField(
107  dtype=float,
108  default=[45000.0],
109  doc="Peak flux level (in DN) of simulated 'astronomical sources'.",
110  )
111  sourceAmp = pexConfig.ListField(
112  dtype=int,
113  default=[0],
114  doc="Amplifier to place simulated 'astronomical sources'.",
115  )
116  sourceX = pexConfig.ListField(
117  dtype=float,
118  default=[50.0],
119  doc="Peak position (in amplifier coordinates) of simulated 'astronomical sources'.",
120  )
121  sourceY = pexConfig.ListField(
122  dtype=float,
123  default=[25.0],
124  doc="Peak position (in amplifier coordinates) of simulated 'astronomical sources'.",
125  )
126  overscanScale = pexConfig.Field(
127  dtype=float,
128  default=100.0,
129  doc="Amplitude (in DN) of the ramp function to add to overscan data.",
130  )
131  biasLevel = pexConfig.Field(
132  dtype=float,
133  default=8000.0,
134  doc="Background contribution to be generated from the bias offset in DN.",
135  )
136  darkRate = pexConfig.Field(
137  dtype=float,
138  default=5.0,
139  doc="Background level contribution (in e-/s) to be generated from dark current.",
140  )
141  darkTime = pexConfig.Field(
142  dtype=float,
143  default=5.0,
144  doc="Exposure time for the dark current contribution.",
145  )
146  flatDrop = pexConfig.Field(
147  dtype=float,
148  default=0.1,
149  doc="Fractional flux drop due to flat from center to edge of detector along x-axis.",
150  )
151  fringeScale = pexConfig.ListField(
152  dtype=float,
153  default=[200.0],
154  doc="Peak fluxes for the components of the fringe ripple in DN.",
155  )
156  fringeX0 = pexConfig.ListField(
157  dtype=float,
158  default=[-100],
159  doc="Center position for the fringe ripples.",
160  )
161  fringeY0 = pexConfig.ListField(
162  dtype=float,
163  default=[-0],
164  doc="Center position for the fringe ripples.",
165  )
166 
167  # Inclusion parameters
168  doAddSky = pexConfig.Field(
169  dtype=bool,
170  default=True,
171  doc="Apply 'sky' signal to output image.",
172  )
173  doAddSource = pexConfig.Field(
174  dtype=bool,
175  default=True,
176  doc="Add simulated source to output image.",
177  )
178  doAddCrosstalk = pexConfig.Field(
179  dtype=bool,
180  default=False,
181  doc="Apply simulated crosstalk to output image. This cannot be corrected by ISR, "
182  "as detector.hasCrosstalk()==False.",
183  )
184  doAddOverscan = pexConfig.Field(
185  dtype=bool,
186  default=True,
187  doc="If untrimmed, add overscan ramp to overscan and data regions.",
188  )
189  doAddBias = pexConfig.Field(
190  dtype=bool,
191  default=True,
192  doc="Add bias signal to data.",
193  )
194  doAddDark = pexConfig.Field(
195  dtype=bool,
196  default=True,
197  doc="Add dark signal to data.",
198  )
199  doAddFlat = pexConfig.Field(
200  dtype=bool,
201  default=True,
202  doc="Add flat signal to data.",
203  )
204  doAddFringe = pexConfig.Field(
205  dtype=bool,
206  default=True,
207  doc="Add fringe signal to data.",
208  )
209 
210  # Datasets to create and return instead of generating an image.
211  doTransmissionCurve = pexConfig.Field(
212  dtype=bool,
213  default=False,
214  doc="Return a simulated transmission curve.",
215  )
216  doDefects = pexConfig.Field(
217  dtype=bool,
218  default=False,
219  doc="Return a simulated defect list.",
220  )
221  doBrighterFatter = pexConfig.Field(
222  dtype=bool,
223  default=False,
224  doc="Return a simulated brighter-fatter kernel.",
225  )
226  doCrosstalkCoeffs = pexConfig.Field(
227  dtype=bool,
228  default=False,
229  doc="Return the matrix of crosstalk coefficients.",
230  )
231  doDataRef = pexConfig.Field(
232  dtype=bool,
233  default=False,
234  doc="Return a simulated gen2 butler dataRef.",
235  )
236  doGenerateImage = pexConfig.Field(
237  dtype=bool,
238  default=False,
239  doc="Return the generated output image if True.",
240  )
241  doGenerateData = pexConfig.Field(
242  dtype=bool,
243  default=False,
244  doc="Return a non-image data structure if True.",
245  )
246  doGenerateAmpDict = pexConfig.Field(
247  dtype=bool,
248  default=False,
249  doc="Return a dict of exposure amplifiers instead of an afwImage.Exposure.",
250  )
251 
252 
253 class IsrMock(pipeBase.Task):
254  """Class to generate consistent mock images for ISR testing.
255 
256  ISR testing currently relies on one-off fake images that do not
257  accurately mimic the full set of detector effects. This class
258  uses the test camera/detector/amplifier structure defined in
259  `lsst.afw.cameraGeom.testUtils` to avoid making the test data
260  dependent on any of the actual obs package formats.
261  """
262  ConfigClass = IsrMockConfig
263  _DefaultName = "isrMock"
264 
265  def __init__(self, **kwargs):
266  super().__init__(**kwargs)
267  self.rngrng = np.random.RandomState(self.config.rngSeed)
268  self.crosstalkCoeffscrosstalkCoeffs = np.array([[0.0, 0.0, 0.0, 0.0, 0.0, -1e-3, 0.0, 0.0],
269  [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
270  [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
271  [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
272  [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
273  [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
274  [1e-2, 0.0, 0.0, 2.2e-2, 0.0, 0.0, 0.0, 0.0],
275  [1e-2, 5e-3, 5e-4, 3e-3, 4e-2, 5e-3, 5e-3, 0.0]])
276 
277  self.bfKernelbfKernel = np.array([[1., 4., 7., 4., 1.],
278  [4., 16., 26., 16., 4.],
279  [7., 26., 41., 26., 7.],
280  [4., 16., 26., 16., 4.],
281  [1., 4., 7., 4., 1.]]) / 273.0
282 
283  def run(self):
284  """Generate a mock ISR product, and return it.
285 
286  Returns
287  -------
288  image : `lsst.afw.image.Exposure`
289  Simulated ISR image with signals added.
290  dataProduct :
291  Simulated ISR data products.
292  None :
293  Returned if no valid configuration was found.
294 
295  Raises
296  ------
297  RuntimeError
298  Raised if both doGenerateImage and doGenerateData are specified.
299  """
300  if self.config.doGenerateImage and self.config.doGenerateData:
301  raise RuntimeError("Only one of doGenerateImage and doGenerateData may be specified.")
302  elif self.config.doGenerateImage:
303  return self.makeImagemakeImage()
304  elif self.config.doGenerateData:
305  return self.makeDatamakeData()
306  else:
307  return None
308 
309  def makeData(self):
310  """Generate simulated ISR data.
311 
312  Currently, only the class defined crosstalk coefficient
313  matrix, brighter-fatter kernel, a constant unity transmission
314  curve, or a simple single-entry defect list can be generated.
315 
316  Returns
317  -------
318  dataProduct :
319  Simulated ISR data product.
320  """
321  if sum(map(bool, [self.config.doBrighterFatter,
322  self.config.doDefects,
323  self.config.doTransmissionCurve,
324  self.config.doCrosstalkCoeffs])) != 1:
325  raise RuntimeError("Only one data product can be generated at a time.")
326  elif self.config.doBrighterFatter is True:
327  return self.makeBfKernelmakeBfKernel()
328  elif self.config.doDefects is True:
329  return self.makeDefectListmakeDefectList()
330  elif self.config.doTransmissionCurve is True:
331  return self.makeTransmissionCurvemakeTransmissionCurve()
332  elif self.config.doCrosstalkCoeffs is True:
333  return self.crosstalkCoeffscrosstalkCoeffs
334  else:
335  return None
336 
337  def makeBfKernel(self):
338  """Generate a simple Gaussian brighter-fatter kernel.
339 
340  Returns
341  -------
342  kernel : `numpy.ndarray`
343  Simulated brighter-fatter kernel.
344  """
345  return self.bfKernelbfKernel
346 
347  def makeDefectList(self):
348  """Generate a simple single-entry defect list.
349 
350  Returns
351  -------
352  defectList : `lsst.meas.algorithms.Defects`
353  Simulated defect list
354  """
356  lsst.geom.Extent2I(40, 50))])
357 
359  """Generate the simulated crosstalk coefficients.
360 
361  Returns
362  -------
363  coeffs : `numpy.ndarray`
364  Simulated crosstalk coefficients.
365  """
366 
367  return self.crosstalkCoeffscrosstalkCoeffs
368 
370  """Generate a simulated flat transmission curve.
371 
372  Returns
373  -------
374  transmission : `lsst.afw.image.TransmissionCurve`
375  Simulated transmission curve.
376  """
377 
378  return afwImage.TransmissionCurve.makeIdentity()
379 
380  def makeImage(self):
381  """Generate a simulated ISR image.
382 
383  Returns
384  -------
385  exposure : `lsst.afw.image.Exposure` or `dict`
386  Simulated ISR image data.
387 
388  Notes
389  -----
390  This method currently constructs a "raw" data image by:
391  * Generating a simulated sky with noise
392  * Adding a single Gaussian "star"
393  * Adding the fringe signal
394  * Multiplying the frame by the simulated flat
395  * Adding dark current (and noise)
396  * Adding a bias offset (and noise)
397  * Adding an overscan gradient parallel to the pixel y-axis
398  * Simulating crosstalk by adding a scaled version of each
399  amplifier to each other amplifier.
400 
401  The exposure with image data constructed this way is in one of
402  three formats.
403  * A single image, with overscan and prescan regions retained
404  * A single image, with overscan and prescan regions trimmed
405  * A `dict`, containing the amplifer data indexed by the
406  amplifier name.
407 
408  The nonlinearity, CTE, and brighter fatter are currently not
409  implemented.
410 
411  Note that this method generates an image in the reverse
412  direction as the ISR processing, as the output image here has
413  had a series of instrument effects added to an idealized
414  exposure.
415  """
416  exposure = self.getExposuregetExposure()
417 
418  for idx, amp in enumerate(exposure.getDetector()):
419  bbox = None
420  if self.config.isTrimmed is True:
421  bbox = amp.getBBox()
422  else:
423  bbox = amp.getRawDataBBox()
424 
425  ampData = exposure.image[bbox]
426 
427  if self.config.doAddSky is True:
428  self.amplifierAddNoiseamplifierAddNoise(ampData, self.config.skyLevel, np.sqrt(self.config.skyLevel))
429 
430  if self.config.doAddSource is True:
431  for sourceAmp, sourceFlux, sourceX, sourceY in zip(self.config.sourceAmp,
432  self.config.sourceFlux,
433  self.config.sourceX,
434  self.config.sourceY):
435  if idx == sourceAmp:
436  self.amplifierAddSourceamplifierAddSource(ampData, sourceFlux, sourceX, sourceY)
437 
438  if self.config.doAddFringe is True:
439  self.amplifierAddFringeamplifierAddFringe(amp, ampData, np.array(self.config.fringeScale),
440  x0=np.array(self.config.fringeX0),
441  y0=np.array(self.config.fringeY0))
442 
443  if self.config.doAddFlat is True:
444  if ampData.getArray().sum() == 0.0:
445  self.amplifierAddNoiseamplifierAddNoise(ampData, 1.0, 0.0)
446  u0 = exposure.getDimensions().getX()
447  v0 = exposure.getDimensions().getY()
448  self.amplifierMultiplyFlatamplifierMultiplyFlat(amp, ampData, self.config.flatDrop, u0=u0, v0=v0)
449 
450  if self.config.doAddDark is True:
451  self.amplifierAddNoiseamplifierAddNoise(ampData,
452  self.config.darkRate * self.config.darkTime / self.config.gain,
453  np.sqrt(self.config.darkRate
454  * self.config.darkTime / self.config.gain))
455 
456  if self.config.doAddCrosstalk is True:
457  ctCalib = CrosstalkCalib()
458  for idxS, ampS in enumerate(exposure.getDetector()):
459  for idxT, ampT in enumerate(exposure.getDetector()):
460  ampDataT = exposure.image[ampT.getBBox()
461  if self.config.isTrimmed else ampT.getRawDataBBox()]
462  outAmp = ctCalib.extractAmp(exposure.getImage(), ampS, ampT,
463  isTrimmed=self.config.isTrimmed)
464  self.amplifierAddCTamplifierAddCT(outAmp, ampDataT, self.crosstalkCoeffscrosstalkCoeffs[idxT][idxS])
465 
466  for amp in exposure.getDetector():
467  bbox = None
468  if self.config.isTrimmed is True:
469  bbox = amp.getBBox()
470  else:
471  bbox = amp.getRawDataBBox()
472 
473  ampData = exposure.image[bbox]
474 
475  if self.config.doAddBias is True:
476  self.amplifierAddNoiseamplifierAddNoise(ampData, self.config.biasLevel,
477  self.config.readNoise / self.config.gain)
478 
479  if self.config.doAddOverscan is True:
480  oscanBBox = amp.getRawHorizontalOverscanBBox()
481  oscanData = exposure.image[oscanBBox]
482  self.amplifierAddNoiseamplifierAddNoise(oscanData, self.config.biasLevel,
483  self.config.readNoise / self.config.gain)
484 
485  self.amplifierAddYGradientamplifierAddYGradient(ampData, -1.0 * self.config.overscanScale,
486  1.0 * self.config.overscanScale)
487  self.amplifierAddYGradientamplifierAddYGradient(oscanData, -1.0 * self.config.overscanScale,
488  1.0 * self.config.overscanScale)
489 
490  if self.config.doGenerateAmpDict is True:
491  expDict = dict()
492  for amp in exposure.getDetector():
493  expDict[amp.getName()] = exposure
494  return expDict
495  else:
496  return exposure
497 
498  # afw primatives to construct the image structure
499  def getCamera(self):
500  """Construct a test camera object.
501 
502  Returns
503  -------
504  camera : `lsst.afw.cameraGeom.camera`
505  Test camera.
506  """
507  cameraWrapper = afwTestUtils.CameraWrapper(
508  plateScale=self.config.plateScale,
509  radialDistortion=self.config.radialDistortion,
510  isLsstLike=self.config.isLsstLike,
511  )
512  camera = cameraWrapper.camera
513  return camera
514 
515  def getExposure(self):
516  """Construct a test exposure.
517 
518  The test exposure has a simple WCS set, as well as a list of
519  unlikely header keywords that can be removed during ISR
520  processing to exercise that code.
521 
522  Returns
523  -------
524  exposure : `lsst.afw.exposure.Exposure`
525  Construct exposure containing masked image of the
526  appropriate size.
527  """
528  camera = self.getCameragetCamera()
529  detector = camera[self.config.detectorIndex]
530  image = afwUtils.makeImageFromCcd(detector,
531  isTrimmed=self.config.isTrimmed,
532  showAmpGain=False,
533  rcMarkSize=0,
534  binSize=1,
535  imageFactory=afwImage.ImageF)
536 
537  var = afwImage.ImageF(image.getDimensions())
538  mask = afwImage.Mask(image.getDimensions())
539  image.assign(0.0)
540 
541  maskedImage = afwImage.makeMaskedImage(image, mask, var)
542  exposure = afwImage.makeExposure(maskedImage)
543  exposure.setDetector(detector)
544  exposure.setWcs(self.getWcsgetWcs())
545 
546  visitInfo = afwImage.VisitInfo(exposureTime=self.config.expTime, darkTime=self.config.darkTime)
547  exposure.getInfo().setVisitInfo(visitInfo)
548 
549  metadata = exposure.getMetadata()
550  metadata.add("SHEEP", 7.3, "number of sheep on farm")
551  metadata.add("MONKEYS", 155, "monkeys per tree")
552  metadata.add("VAMPIRES", 4, "How scary are vampires.")
553 
554  ccd = exposure.getDetector()
555  newCcd = ccd.rebuild()
556  newCcd.clear()
557  for amp in ccd:
558  newAmp = amp.rebuild()
559  newAmp.setLinearityCoeffs((0., 1., 0., 0.))
560  newAmp.setLinearityType("Polynomial")
561  newAmp.setGain(self.config.gain)
562  newAmp.setSuspectLevel(25000.0)
563  newAmp.setSaturation(32000.0)
564  newCcd.append(newAmp)
565  exposure.setDetector(newCcd.finish())
566 
567  exposure.image.array[:] = np.zeros(exposure.getImage().getDimensions()).transpose()
568  exposure.mask.array[:] = np.zeros(exposure.getMask().getDimensions()).transpose()
569  exposure.variance.array[:] = np.zeros(exposure.getVariance().getDimensions()).transpose()
570 
571  return exposure
572 
573  def getWcs(self):
574  """Construct a dummy WCS object.
575 
576  Taken from the deprecated ip_isr/examples/exampleUtils.py.
577 
578  This is not guaranteed, given the distortion and pixel scale
579  listed in the afwTestUtils camera definition.
580 
581  Returns
582  -------
583  wcs : `lsst.afw.geom.SkyWcs`
584  Test WCS transform.
585  """
586  return afwGeom.makeSkyWcs(crpix=lsst.geom.Point2D(0.0, 100.0),
587  crval=lsst.geom.SpherePoint(45.0, 25.0, lsst.geom.degrees),
588  cdMatrix=afwGeom.makeCdMatrix(scale=1.0*lsst.geom.degrees))
589 
590  def localCoordToExpCoord(self, ampData, x, y):
591  """Convert between a local amplifier coordinate and the full
592  exposure coordinate.
593 
594  Parameters
595  ----------
596  ampData : `lsst.afw.image.ImageF`
597  Amplifier image to use for conversions.
598  x : `int`
599  X-coordinate of the point to transform.
600  y : `int`
601  Y-coordinate of the point to transform.
602 
603  Returns
604  -------
605  u : `int`
606  Transformed x-coordinate.
607  v : `int`
608  Transformed y-coordinate.
609 
610  Notes
611  -----
612  The output is transposed intentionally here, to match the
613  internal transpose between numpy and afw.image coordinates.
614  """
615  u = x + ampData.getBBox().getBeginX()
616  v = y + ampData.getBBox().getBeginY()
617 
618  return (v, u)
619 
620  # Simple data values.
621  def amplifierAddNoise(self, ampData, mean, sigma):
622  """Add Gaussian noise to an amplifier's image data.
623 
624  This method operates in the amplifier coordinate frame.
625 
626  Parameters
627  ----------
628  ampData : `lsst.afw.image.ImageF`
629  Amplifier image to operate on.
630  mean : `float`
631  Mean value of the Gaussian noise.
632  sigma : `float`
633  Sigma of the Gaussian noise.
634  """
635  ampArr = ampData.array
636  ampArr[:] = ampArr[:] + self.rngrng.normal(mean, sigma,
637  size=ampData.getDimensions()).transpose()
638 
639  def amplifierAddYGradient(self, ampData, start, end):
640  """Add a y-axis linear gradient to an amplifier's image data.
641 
642  This method operates in the amplifier coordinate frame.
643 
644  Parameters
645  ----------
646  ampData : `lsst.afw.image.ImageF`
647  Amplifier image to operate on.
648  start : `float`
649  Start value of the gradient (at y=0).
650  end : `float`
651  End value of the gradient (at y=ymax).
652  """
653  nPixY = ampData.getDimensions().getY()
654  ampArr = ampData.array
655  ampArr[:] = ampArr[:] + (np.interp(range(nPixY), (0, nPixY - 1), (start, end)).reshape(nPixY, 1)
656  + np.zeros(ampData.getDimensions()).transpose())
657 
658  def amplifierAddSource(self, ampData, scale, x0, y0):
659  """Add a single Gaussian source to an amplifier.
660 
661  This method operates in the amplifier coordinate frame.
662 
663  Parameters
664  ----------
665  ampData : `lsst.afw.image.ImageF`
666  Amplifier image to operate on.
667  scale : `float`
668  Peak flux of the source to add.
669  x0 : `float`
670  X-coordinate of the source peak.
671  y0 : `float`
672  Y-coordinate of the source peak.
673  """
674  for x in range(0, ampData.getDimensions().getX()):
675  for y in range(0, ampData.getDimensions().getY()):
676  ampData.array[y][x] = (ampData.array[y][x]
677  + scale * np.exp(-0.5 * ((x - x0)**2 + (y - y0)**2) / 3.0**2))
678 
679  def amplifierAddCT(self, ampDataSource, ampDataTarget, scale):
680  """Add a scaled copy of an amplifier to another, simulating crosstalk.
681 
682  This method operates in the amplifier coordinate frame.
683 
684  Parameters
685  ----------
686  ampDataSource : `lsst.afw.image.ImageF`
687  Amplifier image to add scaled copy from.
688  ampDataTarget : `lsst.afw.image.ImageF`
689  Amplifier image to add scaled copy to.
690  scale : `float`
691  Flux scale of the copy to add to the target.
692 
693  Notes
694  -----
695  This simulates simple crosstalk between amplifiers.
696  """
697  ampDataTarget.array[:] = (ampDataTarget.array[:]
698  + scale * ampDataSource.array[:])
699 
700  # Functional form data values.
701  def amplifierAddFringe(self, amp, ampData, scale, x0=100, y0=0):
702  """Add a fringe-like ripple pattern to an amplifier's image data.
703 
704  Parameters
705  ----------
706  amp : `~lsst.afw.ampInfo.AmpInfoRecord`
707  Amplifier to operate on. Needed for amp<->exp coordinate transforms.
708  ampData : `lsst.afw.image.ImageF`
709  Amplifier image to operate on.
710  scale : `numpy.array` or `float`
711  Peak intensity scaling for the ripple.
712  x0 : `numpy.array` or `float`, optional
713  Fringe center
714  y0 : `numpy.array` or `float`, optional
715  Fringe center
716 
717  Notes
718  -----
719  This uses an offset sinc function to generate a ripple
720  pattern. True fringes have much finer structure, but this
721  pattern should be visually identifiable. The (x, y)
722  coordinates are in the frame of the amplifier, and (u, v) in
723  the frame of the full trimmed image.
724  """
725  for x in range(0, ampData.getDimensions().getX()):
726  for y in range(0, ampData.getDimensions().getY()):
727  (u, v) = self.localCoordToExpCoordlocalCoordToExpCoord(amp, x, y)
728  ampData.getArray()[y][x] = np.sum((ampData.getArray()[y][x]
729  + scale * np.sinc(((u - x0) / 50)**2
730  + ((v - y0) / 50)**2)))
731 
732  def amplifierMultiplyFlat(self, amp, ampData, fracDrop, u0=100.0, v0=100.0):
733  """Multiply an amplifier's image data by a flat-like pattern.
734 
735  Parameters
736  ----------
737  amp : `lsst.afw.ampInfo.AmpInfoRecord`
738  Amplifier to operate on. Needed for amp<->exp coordinate transforms.
739  ampData : `lsst.afw.image.ImageF`
740  Amplifier image to operate on.
741  fracDrop : `float`
742  Fractional drop from center to edge of detector along x-axis.
743  u0 : `float`
744  Peak location in detector coordinates.
745  v0 : `float`
746  Peak location in detector coordinates.
747 
748  Notes
749  -----
750  This uses a 2-d Gaussian to simulate an illumination pattern
751  that falls off towards the edge of the detector. The (x, y)
752  coordinates are in the frame of the amplifier, and (u, v) in
753  the frame of the full trimmed image.
754  """
755  if fracDrop >= 1.0:
756  raise RuntimeError("Flat fractional drop cannot be greater than 1.0")
757 
758  sigma = u0 / np.sqrt(-2.0 * np.log(fracDrop))
759 
760  for x in range(0, ampData.getDimensions().getX()):
761  for y in range(0, ampData.getDimensions().getY()):
762  (u, v) = self.localCoordToExpCoordlocalCoordToExpCoord(amp, x, y)
763  f = np.exp(-0.5 * ((u - u0)**2 + (v - v0)**2) / sigma**2)
764  ampData.array[y][x] = (ampData.array[y][x] * f)
765 
766 
768  """Generate a raw exposure suitable for ISR.
769  """
770  def __init__(self, **kwargs):
771  super().__init__(**kwargs)
772  self.config.isTrimmed = False
773  self.config.doGenerateImage = True
774  self.config.doGenerateAmpDict = False
775  self.config.doAddOverscan = True
776  self.config.doAddSky = True
777  self.config.doAddSource = True
778  self.config.doAddCrosstalk = False
779  self.config.doAddBias = True
780  self.config.doAddDark = True
781 
782 
784  """Generate a trimmed raw exposure.
785  """
786  def __init__(self, **kwargs):
787  super().__init__(**kwargs)
788  self.config.isTrimmed = True
789  self.config.doAddOverscan = False
790 
791 
793  """Generate a trimmed raw exposure.
794  """
795  def __init__(self, **kwargs):
796  super().__init__(**kwargs)
797  self.config.isTrimmed = True
798  self.config.doGenerateImage = True
799  self.config.doAddOverscan = False
800  self.config.doAddSky = True
801  self.config.doAddSource = True
802  self.config.doAddCrosstalk = False
803 
804  self.config.doAddBias = False
805  self.config.doAddDark = False
806  self.config.doAddFlat = False
807  self.config.doAddFringe = True
808 
809  self.config.biasLevel = 0.0
810  self.config.readNoise = 10.0
811 
812 
814  """Generate a raw exposure dict suitable for ISR.
815  """
816  def __init__(self, **kwargs):
817  super().__init__(**kwargs)
818  self.config.doGenerateAmpDict = True
819 
820 
822  """Parent class for those that make master calibrations.
823  """
824  def __init__(self, **kwargs):
825  super().__init__(**kwargs)
826  self.config.isTrimmed = True
827  self.config.doGenerateImage = True
828  self.config.doAddOverscan = False
829  self.config.doAddSky = False
830  self.config.doAddSource = False
831  self.config.doAddCrosstalk = False
832 
833  self.config.doAddBias = False
834  self.config.doAddDark = False
835  self.config.doAddFlat = False
836  self.config.doAddFringe = False
837 
838 
840  """Simulated master bias calibration.
841  """
842  def __init__(self, **kwargs):
843  super().__init__(**kwargs)
844  self.config.doAddBias = True
845  self.config.readNoise = 10.0
846 
847 
849  """Simulated master dark calibration.
850  """
851  def __init__(self, **kwargs):
852  super().__init__(**kwargs)
853  self.config.doAddDark = True
854  self.config.darkTime = 1.0
855 
856 
858  """Simulated master flat calibration.
859  """
860  def __init__(self, **kwargs):
861  super().__init__(**kwargs)
862  self.config.doAddFlat = True
863 
864 
866  """Simulated master fringe calibration.
867  """
868  def __init__(self, **kwargs):
869  super().__init__(**kwargs)
870  self.config.doAddFringe = True
871 
872 
874  """Simulated untrimmed master fringe calibration.
875  """
876  def __init__(self, **kwargs):
877  super().__init__(**kwargs)
878  self.config.isTrimmed = False
879 
880 
882  """Simulated brighter-fatter kernel.
883  """
884  def __init__(self, **kwargs):
885  super().__init__(**kwargs)
886  self.config.doGenerateImage = False
887  self.config.doGenerateData = True
888  self.config.doBrighterFatter = True
889  self.config.doDefects = False
890  self.config.doCrosstalkCoeffs = False
891  self.config.doTransmissionCurve = False
892 
893 
895  """Simulated defect list.
896  """
897  def __init__(self, **kwargs):
898  super().__init__(**kwargs)
899  self.config.doGenerateImage = False
900  self.config.doGenerateData = True
901  self.config.doBrighterFatter = False
902  self.config.doDefects = True
903  self.config.doCrosstalkCoeffs = False
904  self.config.doTransmissionCurve = False
905 
906 
908  """Simulated crosstalk coefficient matrix.
909  """
910  def __init__(self, **kwargs):
911  super().__init__(**kwargs)
912  self.config.doGenerateImage = False
913  self.config.doGenerateData = True
914  self.config.doBrighterFatter = False
915  self.config.doDefects = False
916  self.config.doCrosstalkCoeffs = True
917  self.config.doTransmissionCurve = False
918 
919 
921  """Simulated transmission curve.
922  """
923  def __init__(self, **kwargs):
924  super().__init__(**kwargs)
925  self.config.doGenerateImage = False
926  self.config.doGenerateData = True
927  self.config.doBrighterFatter = False
928  self.config.doDefects = False
929  self.config.doCrosstalkCoeffs = False
930  self.config.doTransmissionCurve = True
931 
932 
934  """Simulated gen2 butler data ref.
935 
936  Currently only supports get and put operations, which are most
937  likely to be called for data in ISR processing.
938 
939  """
940  dataId = "isrMock Fake Data"
941  darkval = 2. # e-/sec
942  oscan = 250. # DN
943  gradient = .10
944  exptime = 15.0 # seconds
945  darkexptime = 15.0 # seconds
946 
947  def __init__(self, **kwargs):
948  if 'config' in kwargs.keys():
949  self.configconfig = kwargs['config']
950  else:
951  self.configconfig = None
952 
953  def expectImage(self):
954  if self.configconfig is None:
955  self.configconfig = IsrMockConfig()
956  self.configconfig.doGenerateImage = True
957  self.configconfig.doGenerateData = False
958 
959  def expectData(self):
960  if self.configconfig is None:
961  self.configconfig = IsrMockConfig()
962  self.configconfig.doGenerateImage = False
963  self.configconfig.doGenerateData = True
964 
965  def get(self, dataType, **kwargs):
966  """Return an appropriate data product.
967 
968  Parameters
969  ----------
970  dataType : `str`
971  Type of data product to return.
972 
973  Returns
974  -------
975  mock : IsrMock.run() result
976  The output product.
977  """
978  if "_filename" in dataType:
979  self.expectDataexpectData()
980  return tempfile.mktemp(), "mock"
981  elif 'transmission_' in dataType:
982  self.expectDataexpectData()
983  return TransmissionMock(config=self.configconfig).run()
984  elif dataType == 'ccdExposureId':
985  self.expectDataexpectData()
986  return 20090913
987  elif dataType == 'camera':
988  self.expectDataexpectData()
989  return IsrMock(config=self.configconfig).getCamera()
990  elif dataType == 'raw':
991  self.expectImageexpectImage()
992  return RawMock(config=self.configconfig).run()
993  elif dataType == 'bias':
994  self.expectImageexpectImage()
995  return BiasMock(config=self.configconfig).run()
996  elif dataType == 'dark':
997  self.expectImageexpectImage()
998  return DarkMock(config=self.configconfig).run()
999  elif dataType == 'flat':
1000  self.expectImageexpectImage()
1001  return FlatMock(config=self.configconfig).run()
1002  elif dataType == 'fringe':
1003  self.expectImageexpectImage()
1004  return FringeMock(config=self.configconfig).run()
1005  elif dataType == 'defects':
1006  self.expectDataexpectData()
1007  return DefectMock(config=self.configconfig).run()
1008  elif dataType == 'bfKernel':
1009  self.expectDataexpectData()
1010  return BfKernelMock(config=self.configconfig).run()
1011  elif dataType == 'linearizer':
1012  return None
1013  elif dataType == 'crosstalkSources':
1014  return None
1015  else:
1016  raise RuntimeError("ISR DataRefMock cannot return %s.", dataType)
1017 
1018  def put(self, exposure, filename):
1019  """Write an exposure to a FITS file.
1020 
1021  Parameters
1022  ----------
1023  exposure : `lsst.afw.image.Exposure`
1024  Image data to write out.
1025  filename : `str`
1026  Base name of the output file.
1027  """
1028  exposure.writeFits(filename+".fits")
1029 
1030 
1032  """Simulated gen2 butler data ref.
1033 
1034  Currently only supports get and put operations, which are most
1035  likely to be called for data in ISR processing.
1036 
1037  """
1038  dataId = "isrMock Fake Data"
1039  darkval = 2. # e-/sec
1040  oscan = 250. # DN
1041  gradient = .10
1042  exptime = 15 # seconds
1043  darkexptime = 40. # seconds
1044 
1045  def __init__(self, **kwargs):
1046  if 'config' in kwargs.keys():
1047  self.configconfig = kwargs['config']
1048  else:
1049  self.configconfig = IsrMockConfig()
1050  self.configconfig.isTrimmed = True
1051  self.configconfig.doAddFringe = True
1052  self.configconfig.readNoise = 10.0
1053 
1054  def get(self, dataType, **kwargs):
1055  """Return an appropriate data product.
1056 
1057  Parameters
1058  ----------
1059  dataType : `str`
1060  Type of data product to return.
1061 
1062  Returns
1063  -------
1064  mock : IsrMock.run() result
1065  The output product.
1066  """
1067  if "_filename" in dataType:
1068  return tempfile.mktemp(), "mock"
1069  elif 'transmission_' in dataType:
1070  return TransmissionMock(config=self.configconfig).run()
1071  elif dataType == 'ccdExposureId':
1072  return 20090913
1073  elif dataType == 'camera':
1074  return IsrMock(config=self.configconfig).getCamera()
1075  elif dataType == 'raw':
1076  return CalibratedRawMock(config=self.configconfig).run()
1077  elif dataType == 'bias':
1078  return BiasMock(config=self.configconfig).run()
1079  elif dataType == 'dark':
1080  return DarkMock(config=self.configconfig).run()
1081  elif dataType == 'flat':
1082  return FlatMock(config=self.configconfig).run()
1083  elif dataType == 'fringe':
1084  fringes = []
1085  configCopy = copy.deepcopy(self.configconfig)
1086  for scale, x, y in zip(self.configconfig.fringeScale, self.configconfig.fringeX0, self.configconfig.fringeY0):
1087  configCopy.fringeScale = [1.0]
1088  configCopy.fringeX0 = [x]
1089  configCopy.fringeY0 = [y]
1090  fringes.append(FringeMock(config=configCopy).run())
1091  return fringes
1092  elif dataType == 'defects':
1093  return DefectMock(config=self.configconfig).run()
1094  elif dataType == 'bfKernel':
1095  return BfKernelMock(config=self.configconfig).run()
1096  elif dataType == 'linearizer':
1097  return None
1098  elif dataType == 'crosstalkSources':
1099  return None
1100  else:
1101  return None
1102 
1103  def put(self, exposure, filename):
1104  """Write an exposure to a FITS file.
1105 
1106  Parameters
1107  ----------
1108  exposure : `lsst.afw.image.Exposure`
1109  Image data to write out.
1110  filename : `str`
1111  Base name of the output file.
1112  """
1113  exposure.writeFits(filename+".fits")
Represent a 2-dimensional array of bitmask pixels.
Definition: Mask.h:77
Information about a single exposure of an imaging camera.
Definition: VisitInfo.h:68
An integer coordinate rectangle.
Definition: Box.h:55
Point in an unspecified spherical coordinate system.
Definition: SpherePoint.h:57
def __init__(self, **kwargs)
Definition: isrMock.py:884
def __init__(self, **kwargs)
Definition: isrMock.py:842
def __init__(self, **kwargs)
Definition: isrMock.py:851
def put(self, exposure, filename)
Definition: isrMock.py:1018
def get(self, dataType, **kwargs)
Definition: isrMock.py:965
def __init__(self, **kwargs)
Definition: isrMock.py:947
def __init__(self, **kwargs)
Definition: isrMock.py:897
def __init__(self, **kwargs)
Definition: isrMock.py:860
def put(self, exposure, filename)
Definition: isrMock.py:1103
def get(self, dataType, **kwargs)
Definition: isrMock.py:1054
def __init__(self, **kwargs)
Definition: isrMock.py:868
def amplifierAddSource(self, ampData, scale, x0, y0)
Definition: isrMock.py:658
def __init__(self, **kwargs)
Definition: isrMock.py:265
def amplifierAddNoise(self, ampData, mean, sigma)
Definition: isrMock.py:621
def amplifierAddYGradient(self, ampData, start, end)
Definition: isrMock.py:639
def amplifierAddCT(self, ampDataSource, ampDataTarget, scale)
Definition: isrMock.py:679
def amplifierMultiplyFlat(self, amp, ampData, fracDrop, u0=100.0, v0=100.0)
Definition: isrMock.py:732
def amplifierAddFringe(self, amp, ampData, scale, x0=100, y0=0)
Definition: isrMock.py:701
def localCoordToExpCoord(self, ampData, x, y)
Definition: isrMock.py:590
def __init__(self, **kwargs)
Definition: isrMock.py:824
def __init__(self, **kwargs)
Definition: isrMock.py:816
def __init__(self, **kwargs)
Definition: isrMock.py:770
def __init__(self, **kwargs)
Definition: isrMock.py:786
std::shared_ptr< SkyWcs > makeSkyWcs(daf::base::PropertySet &metadata, bool strip=false)
Construct a SkyWcs from FITS keywords.
Definition: SkyWcs.cc:521
Eigen::Matrix2d makeCdMatrix(lsst::geom::Angle const &scale, lsst::geom::Angle const &orientation=0 *lsst::geom::degrees, bool flipX=false)
Make a WCS CD matrix.
Definition: SkyWcs.cc:133
Backwards-compatibility support for depersisting the old Calib (FluxMag0/FluxMag0Err) objects.
std::shared_ptr< Exposure< ImagePixelT, MaskPixelT, VariancePixelT > > makeExposure(MaskedImage< ImagePixelT, MaskPixelT, VariancePixelT > &mimage, std::shared_ptr< geom::SkyWcs const > wcs=std::shared_ptr< geom::SkyWcs const >())
A function to return an Exposure of the correct type (cf.
Definition: Exposure.h:462
MaskedImage< ImagePixelT, MaskPixelT, VariancePixelT > * makeMaskedImage(typename std::shared_ptr< Image< ImagePixelT >> image, typename std::shared_ptr< Mask< MaskPixelT >> mask=Mask< MaskPixelT >(), typename std::shared_ptr< Image< VariancePixelT >> variance=Image< VariancePixelT >())
A function to return a MaskedImage of the correct type (cf.
Definition: MaskedImage.h:1240
def run(self, coaddExposures, bbox, wcs)
Definition: getTemplate.py:603