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