LSST Applications g00274db5b6+edbf708997,g00d0e8bbd7+edbf708997,g199a45376c+5137f08352,g1fd858c14a+1d4b6db739,g262e1987ae+f4d9505c4f,g29ae962dfc+7156fb1a53,g2cef7863aa+73c82f25e4,g35bb328faa+edbf708997,g3e17d7035e+5b3adc59f5,g3fd5ace14f+852fa6fbcb,g47891489e3+6dc8069a4c,g53246c7159+edbf708997,g64539dfbff+9f17e571f4,g67b6fd64d1+6dc8069a4c,g74acd417e5+ae494d68d9,g786e29fd12+af89c03590,g7ae74a0b1c+a25e60b391,g7aefaa3e3d+536efcc10a,g7cc15d900a+d121454f8d,g87389fa792+a4172ec7da,g89139ef638+6dc8069a4c,g8d7436a09f+28c28d8d6d,g8ea07a8fe4+db21c37724,g92c671f44c+9f17e571f4,g98df359435+b2e6376b13,g99af87f6a8+b0f4ad7b8d,gac66b60396+966efe6077,gb88ae4c679+7dec8f19df,gbaa8f7a6c5+38b34f4976,gbf99507273+edbf708997,gc24b5d6ed1+9f17e571f4,gca7fc764a6+6dc8069a4c,gcc769fe2a4+97d0256649,gd7ef33dd92+6dc8069a4c,gdab6d2f7ff+ae494d68d9,gdbb4c4dda9+9f17e571f4,ge410e46f29+6dc8069a4c,geaed405ab2+e194be0d2b,w.2025.47
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", "ElectrostaticBfMock", "DefectMock", "CrosstalkCoeffMock",
26 "TransmissionMock", "MockDataContainer", "MockFringeContainer"]
27
28import copy
29import numpy as np
30import tempfile
31import astropy.time
32
33from datetime import datetime, timezone
34
35import lsst.geom
36import lsst.afw.geom as afwGeom
37import lsst.afw.image as afwImage
38from lsstDebug import getDebugFrame
39
40import lsst.afw.cameraGeom.utils as afwUtils
41import lsst.afw.cameraGeom.testUtils as afwTestUtils
42from lsst.afw.cameraGeom import ReadoutCorner
43import lsst.pex.config as pexConfig
44import lsst.pipe.base as pipeBase
45from .crosstalk import CrosstalkCalib
46from .defects import Defects
47
48
49class IsrMockConfig(pexConfig.Config):
50 """Configuration parameters for isrMock.
51
52 These parameters produce generic fixed position signals from
53 various sources, and combine them in a way that matches how those
54 signals are combined to create real data. The camera used is the
55 test camera defined by the afwUtils code.
56 """
57 # Detector parameters. "Exposure" parameters.
58 isLsstLike = pexConfig.Field(
59 dtype=bool,
60 default=False,
61 doc="If True, products have one raw image per amplifier, otherwise, one raw image per detector.",
62 )
63 plateScale = pexConfig.Field(
64 dtype=float,
65 default=20.0,
66 doc="Plate scale used in constructing mock camera.",
67 )
68 radialDistortion = pexConfig.Field(
69 dtype=float,
70 default=0.925,
71 doc="Radial distortion term used in constructing mock camera.",
72 )
73 isTrimmed = pexConfig.Field(
74 dtype=bool,
75 default=True,
76 doc="If True, amplifiers have been trimmed and mosaicked to remove regions outside the data BBox.",
77 )
78 detectorIndex = pexConfig.Field(
79 dtype=int,
80 default=20,
81 doc="Index for the detector to use. The default value uses a standard 2x4 array of amps.",
82 )
83 rngSeed = pexConfig.Field(
84 dtype=int,
85 default=20000913,
86 doc="Seed for random number generator used to add noise.",
87 )
88 # TODO: DM-18345 Check that mocks scale correctly when gain != 1.0
89 gain = pexConfig.Field(
90 dtype=float,
91 default=1.0,
92 doc="Gain for simulated data in electron/adu.",
93 )
94 readNoise = pexConfig.Field(
95 dtype=float,
96 default=5.0,
97 doc="Read noise of the detector in electron.",
98 )
99 expTime = pexConfig.Field(
100 dtype=float,
101 default=5.0,
102 doc="Exposure time for simulated data.",
103 )
104
105 # Signal parameters
106 skyLevel = pexConfig.Field(
107 dtype=float,
108 default=1000.0,
109 doc="Background contribution to be generated from 'the sky' in "
110 "adu (IsrTask) or electron (IsrTaskLSST).",
111 )
112 sourceFlux = pexConfig.ListField(
113 dtype=float,
114 default=[45000.0],
115 doc="Peak flux level of simulated 'astronomical sources' in "
116 "adu (IsrTask) or electron (IsrTaskLSST).",
117 )
118 sourceAmp = pexConfig.ListField(
119 dtype=int,
120 default=[0],
121 doc="Amplifier to place simulated 'astronomical sources'.",
122 )
123 sourceX = pexConfig.ListField(
124 dtype=float,
125 default=[50.0],
126 doc="Peak position (in amplifier coordinates) of simulated 'astronomical sources'.",
127 )
128 sourceY = pexConfig.ListField(
129 dtype=float,
130 default=[25.0],
131 doc="Peak position (in amplifier coordinates) of simulated 'astronomical sources'.",
132 )
133 overscanScale = pexConfig.Field(
134 dtype=float,
135 default=100.0,
136 doc="Amplitude of the ramp function to add to overscan data in "
137 "adu (IsrTask) or electron (IsrTaskLSST)",
138 )
139 biasLevel = pexConfig.Field(
140 dtype=float,
141 default=8000.0,
142 doc="Background contribution to be generated from the bias offset in adu.",
143 )
144 darkRate = pexConfig.Field(
145 dtype=float,
146 default=5.0,
147 doc="Background level contribution (in electron/s) to be generated from dark current.",
148 )
149 darkTime = pexConfig.Field(
150 dtype=float,
151 default=5.0,
152 doc="Exposure time for the dark current contribution.",
153 )
154 flatDrop = pexConfig.Field(
155 dtype=float,
156 default=0.1,
157 doc="Fractional flux drop due to flat from center to edge of detector along x-axis.",
158 )
159 fringeScale = pexConfig.ListField(
160 dtype=float,
161 default=[200.0],
162 doc="Peak fluxes for the components of the fringe ripple in "
163 "adu (IsrTask) or electron (IsrTaskLSST).",
164 )
165 fringeX0 = pexConfig.ListField(
166 dtype=float,
167 default=[-100],
168 doc="Center position for the fringe ripples.",
169 )
170 fringeY0 = pexConfig.ListField(
171 dtype=float,
172 default=[-0],
173 doc="Center position for the fringe ripples.",
174 )
175
176 # Inclusion parameters
177 doAddSky = pexConfig.Field(
178 dtype=bool,
179 default=True,
180 doc="Apply 'sky' signal to output image.",
181 )
182 doAddSource = pexConfig.Field(
183 dtype=bool,
184 default=True,
185 doc="Add simulated source to output image.",
186 )
187 doAddCrosstalk = pexConfig.Field(
188 dtype=bool,
189 default=False,
190 doc="Apply simulated crosstalk to output image. This cannot be corrected by ISR, "
191 "as detector.hasCrosstalk()==False.",
192 )
193 doAddOverscan = pexConfig.Field(
194 dtype=bool,
195 default=True,
196 doc="If untrimmed, add overscan ramp to overscan and data regions.",
197 )
198 doAddBias = pexConfig.Field(
199 dtype=bool,
200 default=True,
201 doc="Add bias signal to data.",
202 )
203 doAddDark = pexConfig.Field(
204 dtype=bool,
205 default=True,
206 doc="Add dark signal to data.",
207 )
208 doAddFlat = pexConfig.Field(
209 dtype=bool,
210 default=True,
211 doc="Add flat signal to data.",
212 )
213 doAddFringe = pexConfig.Field(
214 dtype=bool,
215 default=True,
216 doc="Add fringe signal to data.",
217 )
218
219 # Datasets to create and return instead of generating an image.
220 doTransmissionCurve = pexConfig.Field(
221 dtype=bool,
222 default=False,
223 doc="Return a simulated transmission curve.",
224 )
225 doDefects = pexConfig.Field(
226 dtype=bool,
227 default=False,
228 doc="Return a simulated defect list.",
229 )
230 doBrighterFatter = pexConfig.Field(
231 dtype=bool,
232 default=False,
233 doc="Return a simulated brighter-fatter kernel.",
234 )
235 brighterFatterCalibType = pexConfig.ChoiceField(
236 dtype=str,
237 doc="Type of calibration to generate if doGenerateData=True, otherwise ignored.",
238 allowed={
239 "KERNEL": "No default suspect values; only config overrides will be used.",
240 "ELECTROSTATIC": "Use the default from the camera model (old defaults).",
241 },
242 default="KERNEL",
243 )
244
245 doDeferredCharge = pexConfig.Field(
246 dtype=bool,
247 default=False,
248 doc="Return a simulated deferred charge calibration.",
249 )
250 doCrosstalkCoeffs = pexConfig.Field(
251 dtype=bool,
252 default=False,
253 doc="Return the matrix of crosstalk coefficients.",
254 )
255 doLinearizer = pexConfig.Field(
256 dtype=bool,
257 default=False,
258 doc="Return linearizer dataset.",
259 )
260 doDataRef = pexConfig.Field(
261 dtype=bool,
262 default=False,
263 doc="Return a simulated gen2 butler dataRef.",
264 )
265 doGenerateImage = pexConfig.Field(
266 dtype=bool,
267 default=False,
268 doc="Return the generated output image if True.",
269 )
270 doGenerateData = pexConfig.Field(
271 dtype=bool,
272 default=False,
273 doc="Return a non-image data structure if True.",
274 )
275 doGenerateAmpDict = pexConfig.Field(
276 dtype=bool,
277 default=False,
278 doc="Return a dict of exposure amplifiers instead of an afwImage.Exposure.",
279 )
280
281
282class IsrMock(pipeBase.Task):
283 """Class to generate consistent mock images for ISR testing.
284
285 ISR testing currently relies on one-off fake images that do not
286 accurately mimic the full set of detector effects. This class
287 uses the test camera/detector/amplifier structure defined in
288 `lsst.afw.cameraGeom.testUtils` to avoid making the test data
289 dependent on any of the actual obs package formats.
290 """
291 ConfigClass = IsrMockConfig
292 _DefaultName = "isrMock"
293
294 def __init__(self, **kwargs):
295 super().__init__(**kwargs)
296 self.rng = np.random.RandomState(self.config.rngSeed)
297 # The coefficients have units adu for IsrTask and have
298 # units electron for IsrTaskLSST.
299 self.crosstalkCoeffs = np.array([[0.0, 0.0, 0.0, 0.0, 0.0, -1e-3, 0.0, 0.0],
300 [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
301 [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
302 [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
303 [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
304 [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
305 [1e-2, 0.0, 0.0, 2.2e-2, 0.0, 0.0, 0.0, 0.0],
306 [1e-2, 5e-3, 5e-4, 3e-3, 4e-2, 5e-3, 5e-3, 0.0]])
307 if getDebugFrame(self._display, "mockCrosstalkCoeffs"):
308 self.crosstalkCoeffs = np.array([[0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
309 [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
310 [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
311 [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
312 [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
313 [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
314 [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
315 [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]])
316 # Generic gaussian BF kernel
317 self.bfKernel = np.array([[1., 4., 7., 4., 1.],
318 [4., 16., 26., 16., 4.],
319 [7., 26., 41., 26., 7.],
320 [4., 16., 26., 16., 4.],
321 [1., 4., 7., 4., 1.]]) / 273.0
322
323 self.aN = np.array([[1., 0.5, 0.0125, 0., 0.],
324 [0.5, 0.010, 0., 0., 0.],
325 [0.0125, 0., 0., 0., 0.],
326 [0, 0., 0., 0., 0.]]) * -1e-7
327 self.aE = np.array([[1., 0.5, 0.0125, 0., 0.],
328 [0.5, 0.010, 0., 0., 0.],
329 [0.0125, 0., 0., 0., 0.],
330 [0, 0., 0., 0., 0.]]) * -1e-7
331 self.aS = self.aN
332 self.aW = self.aE
333 self.aVector = (self.aN, self.aS, self.aE, self.aW)
334
335 def run(self):
336 """Generate a mock ISR product, and return it.
337
338 Returns
339 -------
340 image : `lsst.afw.image.Exposure`
341 Simulated ISR image with signals added.
342 dataProduct :
343 Simulated ISR data products.
344 None :
345 Returned if no valid configuration was found.
346
347 Raises
348 ------
349 RuntimeError
350 Raised if both doGenerateImage and doGenerateData are specified.
351 """
352 if self.config.doGenerateImage and self.config.doGenerateData:
353 raise RuntimeError("Only one of doGenerateImage and doGenerateData may be specified.")
354 elif self.config.doGenerateImage:
355 return self.makeImage()
356 elif self.config.doGenerateData:
357 return self.makeData()
358 else:
359 return None
360
361 def makeData(self):
362 """Generate simulated ISR data.
363
364 Currently, only the class defined crosstalk coefficient
365 matrix, brighter-fatter kernel, a constant unity transmission
366 curve, or a simple single-entry defect list can be generated.
367
368 Returns
369 -------
370 dataProduct :
371 Simulated ISR data product.
372 """
373 if sum(map(bool, [self.config.doBrighterFatter,
374 self.config.doDeferredCharge,
375 self.config.doDefects,
376 self.config.doTransmissionCurve,
377 self.config.doCrosstalkCoeffs,
378 self.config.doLinearizer])) != 1:
379 raise RuntimeError("Only one data product can be generated at a time.")
380 elif self.config.doBrighterFatter:
381 if self.config.brighterFatterCalibType == "KERNEL":
382 return self.makeBfKernel()
383 if self.config.brighterFatterCalibType == "ELECTROSTATIC":
384 return self.makeElectrostaticBf()
385 else:
386 raise ValueError("%s is an unknown Brighter-Fatter calibration type." %
387 self.config.brighterFatterCalibType)
388 elif self.config.doDeferredCharge:
389 return self.makeDeferredChargeCalib()
390 elif self.config.doDefects:
391 return self.makeDefectList()
392 elif self.config.doTransmissionCurve:
393 return self.makeTransmissionCurve()
394 elif self.config.doCrosstalkCoeffs:
395 return self.crosstalkCoeffs
396 elif self.config.doLinearizer:
397 return self.makeLinearizer()
398 else:
399 return None
400
401 def makeBfKernel(self):
402 """Generate a simple Gaussian brighter-fatter kernel.
403
404 Returns
405 -------
406 kernel : `numpy.ndarray`
407 Simulated brighter-fatter kernel.
408 """
409 return self.bfKernel
410
412 """Generate a simple Gaussian brighter-fatter kernel.
413
414 Returns
415 -------
416 kernel : `numpy.ndarray`
417 Simulated brighter-fatter kernel.
418 """
419 return self.aVector
420
422 """Generate a CTI calibration.
423 """
424
425 raise NotImplementedError("Mock deferred charge is not implemented for IsrMock.")
426
427 def makeDefectList(self):
428 """Generate a simple single-entry defect list.
429
430 Returns
431 -------
432 defectList : `lsst.meas.algorithms.Defects`
433 Simulated defect list
434 """
436 lsst.geom.Extent2I(40, 50))])
437
439 """Generate the simulated crosstalk coefficients.
440
441 Returns
442 -------
443 coeffs : `numpy.ndarray`
444 Simulated crosstalk coefficients.
445 """
446
447 return self.crosstalkCoeffs
448
450 """Generate a simulated flat transmission curve.
451
452 Returns
453 -------
454 transmission : `lsst.afw.image.TransmissionCurve`
455 Simulated transmission curve.
456 """
457
458 return afwImage.TransmissionCurve.makeIdentity()
459
460 def makeLinearity(self):
461 """Generate a linearity dataset.
462
463 Returns
464 -------
465 linearizer : `lsst.ip.isr.Linearizer`
466 """
467 raise NotImplementedError("Linearizer not implemented for isrMock.")
468
469 def makeImage(self):
470 """Generate a simulated ISR image.
471
472 Returns
473 -------
474 exposure : `lsst.afw.image.Exposure` or `dict`
475 Simulated ISR image data.
476
477 Notes
478 -----
479 This method currently constructs a "raw" data image by:
480
481 * Generating a simulated sky with noise
482 * Adding a single Gaussian "star"
483 * Adding the fringe signal
484 * Multiplying the frame by the simulated flat
485 * Adding dark current (and noise)
486 * Adding a bias offset (and noise)
487 * Adding an overscan gradient parallel to the pixel y-axis
488 * Simulating crosstalk by adding a scaled version of each
489 amplifier to each other amplifier.
490
491 The exposure with image data constructed this way is in one of
492 three formats.
493
494 * A single image, with overscan and prescan regions retained
495 * A single image, with overscan and prescan regions trimmed
496 * A `dict`, containing the amplifer data indexed by the
497 amplifier name.
498
499 The nonlinearity, CTE, and brighter fatter are currently not
500 implemented.
501
502 Note that this method generates an image in the reverse
503 direction as the ISR processing, as the output image here has
504 had a series of instrument effects added to an idealized
505 exposure.
506 """
507 exposure = self.getExposure()
508
509 for idx, amp in enumerate(exposure.getDetector()):
510 bbox = None
511 if self.config.isTrimmed is True:
512 bbox = amp.getBBox()
513 else:
514 bbox = amp.getRawDataBBox()
515
516 ampData = exposure.image[bbox]
517
518 if self.config.doAddSky is True:
519 self.amplifierAddNoise(ampData, self.config.skyLevel, np.sqrt(self.config.skyLevel))
520
521 if self.config.doAddSource is True:
522 for sourceAmp, sourceFlux, sourceX, sourceY in zip(self.config.sourceAmp,
523 self.config.sourceFlux,
524 self.config.sourceX,
525 self.config.sourceY):
526 if idx == sourceAmp:
527 self.amplifierAddSource(ampData, sourceFlux, sourceX, sourceY)
528
529 if self.config.doAddFringe is True:
530 self.amplifierAddFringe(amp, ampData, np.array(self.config.fringeScale),
531 x0=np.array(self.config.fringeX0),
532 y0=np.array(self.config.fringeY0))
533
534 if self.config.doAddFlat is True:
535 if ampData.getArray().sum() == 0.0:
536 self.amplifierAddNoise(ampData, 1.0, 0.0)
537 u0 = exposure.getDimensions().getX()
538 v0 = exposure.getDimensions().getY()
539 self.amplifierMultiplyFlat(amp, ampData, self.config.flatDrop, u0=u0, v0=v0)
540
541 if self.config.doAddDark is True:
542 self.amplifierAddNoise(ampData,
543 self.config.darkRate * self.config.darkTime / self.config.gain,
544 np.sqrt(self.config.darkRate
545 * self.config.darkTime / self.config.gain))
546
547 if self.config.doAddCrosstalk is True:
548 ctCalib = CrosstalkCalib()
549 # We use the regular subtractCrosstalk code but with a negative
550 # sign on the crosstalk coefficients so it adds instead of
551 # subtracts. We only apply the signal plane (ignoreVariance,
552 # subtrahendMasking) with a very large pixel to mask to ensure
553 # no crosstalk mask bits are set.
554 ctCalib.subtractCrosstalk(
555 exposure,
556 crosstalkCoeffs=-1*self.crosstalkCoeffs,
557 doSubtrahendMasking=True,
558 minPixelToMask=np.inf,
559 ignoreVariance=True,
560 fullAmplifier=False,
561 )
562
563 for amp in exposure.getDetector():
564 bbox = None
565 if self.config.isTrimmed is True:
566 bbox = amp.getBBox()
567 else:
568 bbox = amp.getRawDataBBox()
569
570 ampData = exposure.image[bbox]
571
572 if self.config.doAddBias is True:
573 self.amplifierAddNoise(ampData, self.config.biasLevel,
574 self.config.readNoise / self.config.gain)
575
576 if self.config.doAddOverscan is True:
577 oscanBBox = amp.getRawHorizontalOverscanBBox()
578 oscanData = exposure.image[oscanBBox]
579 self.amplifierAddNoise(oscanData, self.config.biasLevel,
580 self.config.readNoise / self.config.gain)
581
582 self.amplifierAddYGradient(ampData, -1.0 * self.config.overscanScale,
583 1.0 * self.config.overscanScale)
584 self.amplifierAddYGradient(oscanData, -1.0 * self.config.overscanScale,
585 1.0 * self.config.overscanScale)
586
587 if self.config.doGenerateAmpDict is True:
588 expDict = dict()
589 for amp in exposure.getDetector():
590 expDict[amp.getName()] = exposure
591 return expDict
592 else:
593 return exposure
594
595 # afw primatives to construct the image structure
596 def getCamera(self, isForAssembly=False):
597 """Construct a test camera object.
598
599 Parameters
600 -------
601 isForAssembly : `bool`
602 If True, construct a camera with "super raw"
603 orientation (all amplifiers have LL readout
604 corner but still contains the necessary flip
605 and offset info needed for assembly. This is
606 needed if isLsstLike is True. If False, return
607 a camera with bboxes flipped and offset to the
608 correct orientation given the readout corner.
609
610 Returns
611 -------
612 camera : `lsst.afw.cameraGeom.camera`
613 Test camera.
614 """
615 cameraWrapper = afwTestUtils.CameraWrapper(
616 plateScale=self.config.plateScale,
617 radialDistortion=self.config.radialDistortion,
618 isLsstLike=self.config.isLsstLike and isForAssembly,
619 )
620 camera = cameraWrapper.camera
621 return camera
622
623 def getExposure(self, isTrimmed=None):
624 """Construct a test exposure.
625
626 The test exposure has a simple WCS set, as well as a list of
627 unlikely header keywords that can be removed during ISR
628 processing to exercise that code.
629
630 Parameters
631 ----------
632 isTrimmed : `bool` or `None`, optional
633 Override the configuration isTrimmed?
634
635 Returns
636 -------
637 exposure : `lsst.afw.exposure.Exposure`
638 Construct exposure containing masked image of the
639 appropriate size.
640 """
641 if isTrimmed is None:
642 _isTrimmed = self.config.isTrimmed
643 else:
644 _isTrimmed = isTrimmed
645
646 camera = self.getCamera(isForAssembly=self.config.isLsstLike)
647 detector = camera[self.config.detectorIndex]
648 image = afwUtils.makeImageFromCcd(
649 detector,
650 isTrimmed=_isTrimmed,
651 showAmpGain=False,
652 rcMarkSize=0,
653 binSize=1,
654 imageFactory=afwImage.ImageF,
655 )
656
657 var = afwImage.ImageF(image.getDimensions())
658 mask = afwImage.Mask(image.getDimensions())
659 image.assign(0.0)
660
661 maskedImage = afwImage.makeMaskedImage(image, mask, var)
662 exposure = afwImage.makeExposure(maskedImage)
663 exposure.setDetector(detector)
664 exposure.setWcs(self.getWcs())
665
666 visitInfo = afwImage.VisitInfo(exposureTime=self.config.expTime, darkTime=self.config.darkTime)
667 exposure.getInfo().setVisitInfo(visitInfo)
668 # Set a dummy ID.
669 exposure.getInfo().setId(12345)
670
671 metadata = exposure.getMetadata()
672 metadata.add("SHEEP", 7.3, "number of sheep on farm")
673 metadata.add("MONKEYS", 155, "monkeys per tree")
674 metadata.add("VAMPIRES", 4, "How scary are vampires.")
675
676 # Add the current time
677 now = datetime.now(timezone.utc)
678 currentMjd = astropy.time.Time(now, format='datetime', scale='utc').mjd
679 metadata.add("MJD", currentMjd, "Modified Julian Date that the file was written")
680
681 ccd = exposure.getDetector()
682 newCcd = ccd.rebuild()
683 newCcd.clear()
684 readoutMap = {
685 'LL': ReadoutCorner.LL,
686 'LR': ReadoutCorner.LR,
687 'UR': ReadoutCorner.UR,
688 'UL': ReadoutCorner.UL,
689 }
690 for amp in ccd:
691 newAmp = amp.rebuild()
692 newAmp.setLinearityCoeffs((0., 1., 0., 0.))
693 newAmp.setLinearityType("Polynomial")
694 newAmp.setGain(self.config.gain)
695 newAmp.setSuspectLevel(25000.0)
696 newAmp.setSaturation(32000.0)
697 readoutCorner = amp.getReadoutCorner().name
698
699 # Apply flips to bbox where needed
700 imageBBox = amp.getRawDataBBox()
701 rawBbox = amp.getRawBBox()
702 parallelOscanBBox = amp.getRawParallelOverscanBBox()
703 serialOscanBBox = amp.getRawSerialOverscanBBox()
704 prescanBBox = amp.getRawPrescanBBox()
705
706 if self.config.isLsstLike:
707 # This follows cameraGeom.testUtils
708 xoffset, yoffset = amp.getRawXYOffset()
709 offext = lsst.geom.Extent2I(xoffset, yoffset)
710 flipx = bool(amp.getRawFlipX())
711 flipy = bool(amp.getRawFlipY())
712 if flipx:
713 xExt = rawBbox.getDimensions().getX()
714 rawBbox.flipLR(xExt)
715 imageBBox.flipLR(xExt)
716 parallelOscanBBox.flipLR(xExt)
717 serialOscanBBox.flipLR(xExt)
718 prescanBBox.flipLR(xExt)
719 if flipy:
720 yExt = rawBbox.getDimensions().getY()
721 rawBbox.flipTB(yExt)
722 imageBBox.flipTB(yExt)
723 parallelOscanBBox.flipTB(yExt)
724 serialOscanBBox.flipTB(yExt)
725 prescanBBox.flipTB(yExt)
726 if not flipx and not flipy:
727 readoutCorner = 'LL'
728 elif flipx and not flipy:
729 readoutCorner = 'LR'
730 elif flipx and flipy:
731 readoutCorner = 'UR'
732 elif not flipx and flipy:
733 readoutCorner = 'UL'
734 rawBbox.shift(offext)
735 imageBBox.shift(offext)
736 parallelOscanBBox.shift(offext)
737 serialOscanBBox.shift(offext)
738 prescanBBox.shift(offext)
739 newAmp.setReadoutCorner(readoutMap[readoutCorner])
740 newAmp.setRawBBox(rawBbox)
741 newAmp.setRawDataBBox(imageBBox)
742 newAmp.setRawParallelOverscanBBox(parallelOscanBBox)
743 newAmp.setRawSerialOverscanBBox(serialOscanBBox)
744 newAmp.setRawPrescanBBox(prescanBBox)
745 newAmp.setRawFlipX(False)
746 newAmp.setRawFlipY(False)
747 no_offset = lsst.geom.Extent2I(0, 0)
748 newAmp.setRawXYOffset(no_offset)
749
750 newCcd.append(newAmp)
751
752 exposure.setDetector(newCcd.finish())
753
754 exposure.image.array[:] = np.zeros(exposure.getImage().getDimensions()).transpose()
755 exposure.mask.array[:] = np.zeros(exposure.getMask().getDimensions()).transpose()
756 exposure.variance.array[:] = np.zeros(exposure.getVariance().getDimensions()).transpose()
757
758 return exposure
759
760 def getWcs(self):
761 """Construct a dummy WCS object.
762
763 Taken from the deprecated ip_isr/examples/exampleUtils.py.
764
765 This is not guaranteed, given the distortion and pixel scale
766 listed in the afwTestUtils camera definition.
767
768 Returns
769 -------
770 wcs : `lsst.afw.geom.SkyWcs`
771 Test WCS transform.
772 """
773 return afwGeom.makeSkyWcs(crpix=lsst.geom.Point2D(0.0, 100.0),
774 crval=lsst.geom.SpherePoint(45.0, 25.0, lsst.geom.degrees),
775 cdMatrix=afwGeom.makeCdMatrix(scale=1.0*lsst.geom.degrees))
776
777 def localCoordToExpCoord(self, ampData, x, y):
778 """Convert between a local amplifier coordinate and the full
779 exposure coordinate.
780
781 Parameters
782 ----------
783 ampData : `lsst.afw.image.ImageF`
784 Amplifier image to use for conversions.
785 x : `int`
786 X-coordinate of the point to transform.
787 y : `int`
788 Y-coordinate of the point to transform.
789
790 Returns
791 -------
792 u : `int`
793 Transformed x-coordinate.
794 v : `int`
795 Transformed y-coordinate.
796
797 Notes
798 -----
799 The output is transposed intentionally here, to match the
800 internal transpose between numpy and afw.image coordinates.
801 """
802 u = x + ampData.getBBox().getBeginX()
803 v = y + ampData.getBBox().getBeginY()
804
805 return (v, u)
806
807 # Simple data values.
808 def amplifierAddNoise(self, ampData, mean, sigma, rng=None):
809 """Add Gaussian noise to an amplifier's image data.
810
811 This method operates in the amplifier coordinate frame.
812
813 Parameters
814 ----------
815 ampData : `lsst.afw.image.ImageF`
816 Amplifier image to operate on.
817 mean : `float`
818 Mean value of the Gaussian noise.
819 sigma : `float`
820 Sigma of the Gaussian noise.
821 rng : `np.random.RandomState`, optional
822 Random state to use instead of self.rng.
823 """
824 if rng is not None:
825 _rng = rng
826 else:
827 _rng = self.rng
828
829 ampArr = ampData.array
830 ampArr[:] = ampArr[:] + _rng.normal(mean, sigma,
831 size=ampData.getDimensions()).transpose()
832
833 def amplifierAddYGradient(self, ampData, start, end):
834 """Add a y-axis linear gradient to an amplifier's image data.
835
836 This method operates in the amplifier coordinate frame.
837
838 Parameters
839 ----------
840 ampData : `lsst.afw.image.ImageF`
841 Amplifier image to operate on.
842 start : `float`
843 Start value of the gradient (at y=0).
844 end : `float`
845 End value of the gradient (at y=ymax).
846 """
847 nPixY = ampData.getDimensions().getY()
848 ampArr = ampData.array
849 ampArr[:] = ampArr[:] + (np.interp(range(nPixY), (0, nPixY - 1), (start, end)).reshape(nPixY, 1)
850 + np.zeros(ampData.getDimensions()).transpose())
851
852 def amplifierAddSource(self, ampData, scale, x0, y0):
853 """Add a single Gaussian source to an amplifier.
854
855 This method operates in the amplifier coordinate frame.
856
857 Parameters
858 ----------
859 ampData : `lsst.afw.image.ImageF`
860 Amplifier image to operate on.
861 scale : `float`
862 Peak flux of the source to add.
863 x0 : `float`
864 X-coordinate of the source peak.
865 y0 : `float`
866 Y-coordinate of the source peak.
867 """
868 for x in range(0, ampData.getDimensions().getX()):
869 for y in range(0, ampData.getDimensions().getY()):
870 ampData.array[y][x] = (ampData.array[y][x]
871 + scale * np.exp(-0.5 * ((x - x0)**2 + (y - y0)**2) / 3.0**2))
872
873 # Functional form data values.
874 def amplifierAddFringe(self, amp, ampData, scale, x0=100, y0=0):
875 """Add a fringe-like ripple pattern to an amplifier's image data.
876
877 Parameters
878 ----------
879 amp : `~lsst.afw.ampInfo.AmpInfoRecord`
880 Amplifier to operate on. Needed for amp<->exp coordinate
881 transforms.
882 ampData : `lsst.afw.image.ImageF`
883 Amplifier image to operate on.
884 scale : `numpy.array` or `float`
885 Peak intensity scaling for the ripple.
886 x0 : `numpy.array` or `float`, optional
887 Fringe center
888 y0 : `numpy.array` or `float`, optional
889 Fringe center
890
891 Notes
892 -----
893 This uses an offset sinc function to generate a ripple
894 pattern. True fringes have much finer structure, but this
895 pattern should be visually identifiable. The (x, y)
896 coordinates are in the frame of the amplifier, and (u, v) in
897 the frame of the full trimmed image.
898 """
899 for x in range(0, ampData.getDimensions().getX()):
900 for y in range(0, ampData.getDimensions().getY()):
901 (u, v) = self.localCoordToExpCoord(amp, x, y)
902 ampData.getArray()[y][x] = np.sum((ampData.getArray()[y][x]
903 + scale * np.sinc(((u - x0) / 50)**2
904 + ((v - y0) / 50)**2)))
905
906 def amplifierMultiplyFlat(self, amp, ampData, fracDrop, u0=100.0, v0=100.0):
907 """Multiply an amplifier's image data by a flat-like pattern.
908
909 Parameters
910 ----------
911 amp : `lsst.afw.ampInfo.AmpInfoRecord`
912 Amplifier to operate on. Needed for amp<->exp coordinate
913 transforms.
914 ampData : `lsst.afw.image.ImageF`
915 Amplifier image to operate on.
916 fracDrop : `float`
917 Fractional drop from center to edge of detector along x-axis.
918 u0 : `float`
919 Peak location in detector coordinates.
920 v0 : `float`
921 Peak location in detector coordinates.
922
923 Notes
924 -----
925 This uses a 2-d Gaussian to simulate an illumination pattern
926 that falls off towards the edge of the detector. The (x, y)
927 coordinates are in the frame of the amplifier, and (u, v) in
928 the frame of the full trimmed image.
929 """
930 if fracDrop >= 1.0:
931 raise RuntimeError("Flat fractional drop cannot be greater than 1.0")
932
933 sigma = u0 / np.sqrt(-2.0 * np.log(fracDrop))
934
935 for x in range(0, ampData.getDimensions().getX()):
936 for y in range(0, ampData.getDimensions().getY()):
937 (u, v) = self.localCoordToExpCoord(amp, x, y)
938 f = np.exp(-0.5 * ((u - u0)**2 + (v - v0)**2) / sigma**2)
939 ampData.array[y][x] = (ampData.array[y][x] * f)
940
941
943 """Generate a raw exposure suitable for ISR.
944 """
945 def __init__(self, **kwargs):
946 super().__init__(**kwargs)
947 self.config.isTrimmed = False
948 self.config.doGenerateImage = True
949 self.config.doGenerateAmpDict = False
950 self.config.doAddOverscan = True
951 self.config.doAddSky = True
952 self.config.doAddSource = True
953 self.config.doAddCrosstalk = False
954 self.config.doAddBias = True
955 self.config.doAddDark = True
956
957
959 """Generate a trimmed raw exposure.
960 """
961 def __init__(self, **kwargs):
962 super().__init__(**kwargs)
963 self.config.isTrimmed = True
964 self.config.doAddOverscan = False
965
966
968 """Generate a trimmed raw exposure.
969 """
970 def __init__(self, **kwargs):
971 super().__init__(**kwargs)
972 self.config.isTrimmed = True
973 self.config.doGenerateImage = True
974 self.config.doAddOverscan = False
975 self.config.doAddSky = True
976 self.config.doAddSource = True
977 self.config.doAddCrosstalk = False
978
979 self.config.doAddBias = False
980 self.config.doAddDark = False
981 self.config.doAddFlat = False
982 self.config.doAddFringe = True
983
984 self.config.biasLevel = 0.0
985 self.config.readNoise = 10.0
986
987
989 """Generate a raw exposure dict suitable for ISR.
990 """
991 def __init__(self, **kwargs):
992 super().__init__(**kwargs)
993 self.config.doGenerateAmpDict = True
994
995
997 """Parent class for those that make master calibrations.
998 """
999 def __init__(self, **kwargs):
1000 super().__init__(**kwargs)
1001 self.config.isTrimmed = True
1002 self.config.doGenerateImage = True
1003 self.config.doAddOverscan = False
1004 self.config.doAddSky = False
1005 self.config.doAddSource = False
1006 self.config.doAddCrosstalk = False
1007
1008 self.config.doAddBias = False
1009 self.config.doAddDark = False
1010 self.config.doAddFlat = False
1011 self.config.doAddFringe = False
1012
1013
1015 """Simulated master bias calibration.
1016 """
1017 def __init__(self, **kwargs):
1018 super().__init__(**kwargs)
1019 self.config.doAddBias = True
1020 self.config.readNoise = 10.0
1021
1022
1024 """Simulated master dark calibration.
1025 """
1026 def __init__(self, **kwargs):
1027 super().__init__(**kwargs)
1028 self.config.doAddDark = True
1029 self.config.darkTime = 1.0
1030
1031
1033 """Simulated master flat calibration.
1034 """
1035 def __init__(self, **kwargs):
1036 super().__init__(**kwargs)
1037 self.config.doAddFlat = True
1038
1039
1041 """Simulated master fringe calibration.
1042 """
1043 def __init__(self, **kwargs):
1044 super().__init__(**kwargs)
1045 self.config.doAddFringe = True
1046
1047
1049 """Simulated untrimmed master fringe calibration.
1050 """
1051 def __init__(self, **kwargs):
1052 super().__init__(**kwargs)
1053 self.config.isTrimmed = False
1054
1055
1057 """Simulated brighter-fatter kernel.
1058 """
1059 def __init__(self, **kwargs):
1060 super().__init__(**kwargs)
1061 self.config.doGenerateImage = False
1062 self.config.doGenerateData = True
1063 self.config.doBrighterFatter = True
1064 self.config.doDefects = False
1065 self.config.doCrosstalkCoeffs = False
1066 self.config.doTransmissionCurve = False
1067
1068
1070 """Simulated brighter-fatter kernel.
1071 """
1072 def __init__(self, **kwargs):
1073 super().__init__(**kwargs)
1074 self.config.doGenerateImage = False
1075 self.config.doGenerateData = True
1076 self.config.doBrighterFatter = True
1077 self.config.doDefects = False
1078 self.config.doCrosstalkCoeffs = False
1079 self.config.doTransmissionCurve = False
1080
1081
1083 """Simulated deferred charge calibration.
1084 """
1085 def __init__(self, **kwargs):
1086 super().__init__(**kwargs)
1087 self.config.doGenerateImage = False
1088 self.config.doGenerateData = True
1089 self.config.doDeferredCharge = True
1090 self.config.doDefects = False
1091 self.config.doCrosstalkCoeffs = False
1092 self.config.doTransmissionCurve = False
1093
1094
1096 """Simulated defect list.
1097 """
1098 def __init__(self, **kwargs):
1099 super().__init__(**kwargs)
1100 self.config.doGenerateImage = False
1101 self.config.doGenerateData = True
1102 self.config.doBrighterFatter = False
1103 self.config.doDefects = True
1104 self.config.doCrosstalkCoeffs = False
1105 self.config.doTransmissionCurve = False
1106
1107
1109 """Simulated crosstalk coefficient matrix.
1110 """
1111 def __init__(self, **kwargs):
1112 super().__init__(**kwargs)
1113 self.config.doGenerateImage = False
1114 self.config.doGenerateData = True
1115 self.config.doBrighterFatter = False
1116 self.config.doDefects = False
1117 self.config.doCrosstalkCoeffs = True
1118 self.config.doTransmissionCurve = False
1119
1120
1122 """Simulated transmission curve.
1123 """
1124 def __init__(self, **kwargs):
1125 super().__init__(**kwargs)
1126 self.config.doGenerateImage = False
1127 self.config.doGenerateData = True
1128 self.config.doBrighterFatter = False
1129 self.config.doDefects = False
1130 self.config.doCrosstalkCoeffs = False
1131 self.config.doTransmissionCurve = True
1132
1133
1134class MockDataContainer(object):
1135 """Container for holding ISR mock objects.
1136 """
1137 dataId = "isrMock Fake Data"
1138 darkval = 2. # electron/sec
1139 oscan = 250. # adu
1140 gradient = .10
1141 exptime = 15.0 # seconds
1142 darkexptime = 15.0 # seconds
1143
1144 def __init__(self, **kwargs):
1145 if 'config' in kwargs.keys():
1146 self.config = kwargs['config']
1147 else:
1148 self.config = None
1149
1150 def expectImage(self):
1151 if self.config is None:
1152 self.config = IsrMockConfig()
1153 self.config.doGenerateImage = True
1154 self.config.doGenerateData = False
1155
1156 def expectData(self):
1157 if self.config is None:
1158 self.config = IsrMockConfig()
1159 self.config.doGenerateImage = False
1160 self.config.doGenerateData = True
1161
1162 def get(self, dataType, **kwargs):
1163 """Return an appropriate data product.
1164
1165 Parameters
1166 ----------
1167 dataType : `str`
1168 Type of data product to return.
1169
1170 Returns
1171 -------
1172 mock : IsrMock.run() result
1173 The output product.
1174 """
1175 if "_filename" in dataType:
1176 self.expectData()
1177 return tempfile.mktemp(), "mock"
1178 elif 'transmission_' in dataType:
1179 self.expectData()
1180 return TransmissionMock(config=self.config).run()
1181 elif dataType == 'ccdExposureId':
1182 self.expectData()
1183 return 20090913
1184 elif dataType == 'camera':
1185 self.expectData()
1186 return IsrMock(config=self.config).getCamera()
1187 elif dataType == 'raw':
1188 self.expectImage()
1189 return RawMock(config=self.config).run()
1190 elif dataType == 'bias':
1191 self.expectImage()
1192 return BiasMock(config=self.config).run()
1193 elif dataType == 'dark':
1194 self.expectImage()
1195 return DarkMock(config=self.config).run()
1196 elif dataType == 'flat':
1197 self.expectImage()
1198 return FlatMock(config=self.config).run()
1199 elif dataType == 'fringe':
1200 self.expectImage()
1201 return FringeMock(config=self.config).run()
1202 elif dataType == 'defects':
1203 self.expectData()
1204 return DefectMock(config=self.config).run()
1205 elif dataType == 'bfKernel':
1206 self.expectData()
1207 return BfKernelMock(config=self.config).run()
1208 elif dataType == 'ebf':
1209 self.expectData()
1210 return ElectrostaticBfMock(config=self.config).run()
1211 elif dataType == 'linearizer':
1212 return None
1213 elif dataType == 'crosstalkSources':
1214 return None
1215 else:
1216 raise RuntimeError("ISR DataRefMock cannot return %s.", dataType)
1217
1218
1220 """Container for mock fringe data.
1221 """
1222 dataId = "isrMock Fake Data"
1223 darkval = 2. # electron/sec
1224 oscan = 250. # adu
1225 gradient = .10
1226 exptime = 15 # seconds
1227 darkexptime = 40. # seconds
1228
1229 def __init__(self, **kwargs):
1230 if 'config' in kwargs.keys():
1231 self.config = kwargs['config']
1232 else:
1233 self.config = IsrMockConfig()
1234 self.config.isTrimmed = True
1235 self.config.doAddFringe = True
1236 self.config.readNoise = 10.0
1237
1238 def get(self, dataType, **kwargs):
1239 """Return an appropriate data product.
1240
1241 Parameters
1242 ----------
1243 dataType : `str`
1244 Type of data product to return.
1245
1246 Returns
1247 -------
1248 mock : IsrMock.run() result
1249 The output product.
1250 """
1251 if "_filename" in dataType:
1252 return tempfile.mktemp(), "mock"
1253 elif 'transmission_' in dataType:
1254 return TransmissionMock(config=self.config).run()
1255 elif dataType == 'ccdExposureId':
1256 return 20090913
1257 elif dataType == 'camera':
1258 return IsrMock(config=self.config).getCamera()
1259 elif dataType == 'raw':
1260 return CalibratedRawMock(config=self.config).run()
1261 elif dataType == 'bias':
1262 return BiasMock(config=self.config).run()
1263 elif dataType == 'dark':
1264 return DarkMock(config=self.config).run()
1265 elif dataType == 'flat':
1266 return FlatMock(config=self.config).run()
1267 elif dataType == 'fringe':
1268 fringes = []
1269 configCopy = copy.deepcopy(self.config)
1270 for scale, x, y in zip(self.config.fringeScale, self.config.fringeX0, self.config.fringeY0):
1271 configCopy.fringeScale = [1.0]
1272 configCopy.fringeX0 = [x]
1273 configCopy.fringeY0 = [y]
1274 fringes.append(FringeMock(config=configCopy).run())
1275 return fringes
1276 elif dataType == 'defects':
1277 return DefectMock(config=self.config).run()
1278 elif dataType == 'bfKernel':
1279 return BfKernelMock(config=self.config).run()
1280 elif dataType == 'ebf':
1281 return ElectrostaticBfMock(config=self.config).run()
1282 elif dataType == 'linearizer':
1283 return None
1284 elif dataType == 'crosstalkSources':
1285 return None
1286 else:
1287 return None
Represent a 2-dimensional array of bitmask pixels.
Definition Mask.h:82
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:777
amplifierAddNoise(self, ampData, mean, sigma, rng=None)
Definition isrMock.py:808
amplifierAddSource(self, ampData, scale, x0, y0)
Definition isrMock.py:852
amplifierMultiplyFlat(self, amp, ampData, fracDrop, u0=100.0, v0=100.0)
Definition isrMock.py:906
amplifierAddYGradient(self, ampData, start, end)
Definition isrMock.py:833
getCamera(self, isForAssembly=False)
Definition isrMock.py:596
getExposure(self, isTrimmed=None)
Definition isrMock.py:623
amplifierAddFringe(self, amp, ampData, scale, x0=100, y0=0)
Definition isrMock.py:874
__init__(self, **kwargs)
Definition isrMock.py:294
get(self, dataType, **kwargs)
Definition isrMock.py:1162
get(self, dataType, **kwargs)
Definition isrMock.py:1238
__init__(self, **kwargs)
Definition isrMock.py:945
std::shared_ptr< SkyWcs > makeSkyWcs(daf::base::PropertySet &metadata, bool strip=false)
Construct a SkyWcs from FITS keywords.
Definition SkyWcs.cc:545
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:150
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:459