Coverage for python / lsst / ip / isr / isrMockLSST.py: 16%
489 statements
« prev ^ index » next coverage.py v7.13.5, created at 2026-04-24 08:28 +0000
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/>.
22__all__ = ["IsrMockLSSTConfig", "IsrMockLSST", "RawMockLSST",
23 "CalibratedRawMockLSST", "ReferenceMockLSST",
24 "BiasMockLSST", "DarkMockLSST", "FlatMockLSST", "FringeMockLSST",
25 "BfKernelMockLSST", "ElectrostaticBfMockLSST", "DefectMockLSST",
26 "CrosstalkCoeffMockLSST", "TransmissionMockLSST"]
28import numpy as np
29import galsim
31import lsst.geom as geom
32import lsst.pex.config as pexConfig
33import lsst.afw.detection as afwDetection
34import lsst.afw.math as afwMath
35from lsst.afw.cameraGeom import ReadoutCorner
36from lsst.daf.base import PropertyList
37from .crosstalk import CrosstalkCalib
38from .isrMock import IsrMockConfig, IsrMock
39from .defects import Defects
40from .assembleCcdTask import AssembleCcdTask
41from .linearize import Linearizer
42from .brighterFatterKernel import BrighterFatterKernel
43from .electrostaticBrighterFatter import ElectrostaticBrighterFatterDistortionMatrix
44from .deferredCharge import (
45 SegmentSimulator,
46 FloatingOutputAmplifier,
47 DeferredChargeCalib
48)
51class IsrMockLSSTConfig(IsrMockConfig):
52 """Configuration parameters for isrMockLSST.
53 """
54 # Detector parameters and "Exposure" parameters,
55 # mostly inherited from IsrMockConfig.
56 isLsstLike = pexConfig.Field(
57 dtype=bool,
58 default=True,
59 doc="If True, products have one raw image per amplifier, otherwise, one raw image per detector.",
60 )
61 calibMode = pexConfig.Field(
62 dtype=bool,
63 default=False,
64 doc="Set to true to produce mock calibration products, e.g. combined bias, dark, flat, etc.",
65 )
66 doAdd2DBias = pexConfig.Field(
67 dtype=bool,
68 default=True,
69 doc="Add 2D bias residual frame to data.",
70 )
71 doAddBrightDefects = pexConfig.Field(
72 dtype=bool,
73 default=True,
74 doc="Add bright defects (bad column) to data.",
75 )
76 brightDefectLevel = pexConfig.Field(
77 dtype=float,
78 default=30000.0,
79 doc="Bright defect level (electron).",
80 )
81 doAddBadParallelOverscanColumn = pexConfig.Field(
82 dtype=bool,
83 default=True,
84 doc="Add a bad column to the parallel overscan.",
85 )
86 badParallelOverscanColumnLevel = pexConfig.Field(
87 dtype=float,
88 default=300000.,
89 doc="Bright parallel overscan column level (electron). Should be above saturation.",
90 )
91 doAddBadParallelOverscanColumnNeighbors = pexConfig.Field(
92 dtype=bool,
93 default=True,
94 doc="Add low-level bad columns next to parallel overscan bad column.",
95 )
96 badParallelOverscanColumnNeighborsLevel = pexConfig.Field(
97 dtype=float,
98 default=50.0,
99 doc="Bright parallel overscan column neighbors level (electron).",
100 )
101 doAddBrighterFatter = pexConfig.Field(
102 dtype=bool,
103 default=False,
104 doc="Add brighter fatter and/or diffusion effects to image.",
105 )
106 doAddDeferredCharge = pexConfig.Field(
107 dtype=bool,
108 default=False,
109 doc="Add serial CTI at the amp level?",
110 )
111 bfStrength = pexConfig.Field(
112 dtype=float,
113 default=2.0,
114 doc="The brighter fatter effect scaling parameter (cannot be zero)."
115 "Nominally = 1, but = 2 is more realistic."
116 )
117 nRecalc = pexConfig.Field(
118 dtype=int,
119 default=10000,
120 doc="Number of electrons to accumulate before recalculating pixel shapes.",
121 )
122 doAddClockInjectedOffset = pexConfig.Field(
123 dtype=bool,
124 default=True,
125 doc="Add clock-injected offset to data (on-chip bias level).",
126 )
127 clockInjectedOffsetLevel = pexConfig.Field(
128 dtype=float,
129 default=8500.0,
130 doc="Clock-injected offset (on-chip bias level), in electron.",
131 )
132 noise2DBias = pexConfig.Field(
133 dtype=float,
134 default=2.0,
135 doc="Noise (in electron) to generate a 2D bias residual frame.",
136 )
137 doAddDarkNoiseOnly = pexConfig.Field(
138 dtype=bool,
139 default=False,
140 doc="Add only dark current noise, for testing consistency.",
141 )
142 doAddParallelOverscanRamp = pexConfig.Field(
143 dtype=bool,
144 default=True,
145 doc="Add overscan ramp to parallel overscan and data regions.",
146 )
147 doAddSerialOverscanRamp = pexConfig.Field(
148 dtype=bool,
149 default=True,
150 doc="Add overscan ramp to serial overscan and data regions.",
151 )
152 doAddHighSignalNonlinearity = pexConfig.Field(
153 dtype=bool,
154 default=True,
155 doc="Add high signal non-linearity to overscan and data regions?",
156 )
157 doAddLowSignalNonlinearity = pexConfig.Field(
158 dtype=bool,
159 default=False,
160 doc="Add low signal non-linearity to overscan and data regions? (Not supported yet.)",
161 )
162 highSignalNonlinearityThreshold = pexConfig.Field(
163 dtype=float,
164 default=40_000.,
165 doc="Threshold (in adu) for the non-linearity to be considered ``high signal``.",
166 )
167 doApplyGain = pexConfig.Field(
168 dtype=bool,
169 default=True,
170 doc="Add gain to data.",
171 )
172 doRoundAdu = pexConfig.Field(
173 dtype=bool,
174 default=True,
175 doc="Round adu values to nearest integer.",
176 )
177 gainDict = pexConfig.DictField(
178 keytype=str,
179 itemtype=float,
180 doc="Dictionary of amp name to gain; any amps not listed will use "
181 "config.gain as the value. Units are electron/adu.",
182 default={
183 "C:0,0": 1.65,
184 "C:0,1": 1.60,
185 "C:0,2": 1.55,
186 "C:0,3": 1.70,
187 "C:1,0": 1.75,
188 "C:1,1": 1.80,
189 "C:1,2": 1.85,
190 "C:1,3": 1.70,
191 },
192 )
193 assembleCcd = pexConfig.ConfigurableField(
194 target=AssembleCcdTask,
195 doc="CCD assembly task; used for defect box conversions.",
196 )
198 def validate(self):
199 super().validate()
201 if self.doAddLowSignalNonlinearity:
202 raise NotImplementedError("Low signal non-linearity is not implemented.")
204 if self.doAddDeferredCharge and self.isTrimmed:
205 raise NotImplementedError("Must be untrimmed for mock serial CTI to"
206 "realistically add charge into overscan regions.")
208 def setDefaults(self):
209 super().setDefaults()
211 self.gain = 1.7 # Default value.
212 self.skyLevel = 1700.0 # electron
213 self.sourceFlux = [50_000.0] # electron
214 self.sourceX = [35.0] # pixel
215 self.sourceY = [37.0] # pixel
216 self.overscanScale = 170.0 # electron
217 self.biasLevel = 20_000.0 # adu
218 self.doAddCrosstalk = True
221class IsrMockLSST(IsrMock):
222 """Class to generate consistent mock images for ISR testing.
223 """
224 ConfigClass = IsrMockLSSTConfig
225 _DefaultName = "isrMockLSST"
227 def __init__(self, **kwargs):
228 super().__init__(**kwargs)
230 # Get kernel derived from imSim generated flats with BFE. The kernel
231 # was used for Ops Rehearsal 3 for LSSTCam-type sensors
232 # See https://rubinobs.atlassian.net/browse/DM-43059 for more details.
233 self.bfKernel = np.array([[4.83499829e-01, 8.10171823e-01, 5.31096720e-01,
234 3.54369868e-02, -8.44782871e-01, -1.64614462e+00,
235 -3.83933101e+00, -5.60243416e+00, -6.51691578e+00,
236 -5.60243416e+00, -3.83933101e+00, -1.64614462e+00,
237 -8.44782871e-01, 3.54369868e-02, 5.31096720e-01,
238 8.10171823e-01, 4.83499829e-01],
239 [1.12382749e+00, 2.22609074e+00, 1.27877807e+00,
240 4.55434098e-01, -1.76842385e+00, -1.90046460e+00,
241 -8.10874526e+00, -1.20534899e+01, -1.48627948e+01,
242 -1.20534899e+01, -8.10874526e+00, -1.90046460e+00,
243 -1.76842385e+00, 4.55434098e-01, 1.27877807e+00,
244 2.22609074e+00, 1.12382749e+00],
245 [1.78571940e+00, 4.38918110e+00, 3.95098587e+00,
246 3.70961649e-01, -3.48151981e+00, -9.61567736e+00,
247 -1.78621172e+01, -2.32278872e+01, -2.31833727e+01,
248 -2.32278872e+01, -1.78621172e+01, -9.61567736e+00,
249 -3.48151981e+00, 3.70961649e-01, 3.95098587e+00,
250 4.38918110e+00, 1.78571940e+00],
251 [1.62986900e+00, 3.67851228e+00, 5.68645252e+00,
252 2.15342566e-01, -8.89937202e+00, -1.44739813e+01,
253 -2.98952660e+01, -4.37420817e+01, -4.83160958e+01,
254 -4.37420817e+01, -2.98952660e+01, -1.44739813e+01,
255 -8.89937202e+00, 2.15342566e-01, 5.68645252e+00,
256 3.67851228e+00, 1.62986900e+00],
257 [1.05524430e+00, 1.71917897e+00, 1.73105590e+00,
258 -2.10088420e+00, -1.15118208e+01, -2.55007598e+01,
259 -4.73056159e+01, -6.97257685e+01, -8.09264433e+01,
260 -6.97257685e+01, -4.73056159e+01, -2.55007598e+01,
261 -1.15118208e+01, -2.10088420e+00, 1.73105590e+00,
262 1.71917897e+00, 1.05524430e+00],
263 [8.71929228e-01, 5.41025574e-01, 9.47560771e-01,
264 -5.75314708e-01, -7.46104027e+00, -4.42314481e+01,
265 -9.54126971e+01, -1.61603201e+02, -2.07520692e+02,
266 -1.61603201e+02, -9.54126971e+01, -4.42314481e+01,
267 -7.46104027e+00, -5.75314708e-01, 9.47560771e-01,
268 5.41025574e-01, 8.71929228e-01],
269 [1.89144704e+00, 3.57543979e+00, -6.91419168e-02,
270 -3.37950835e+00, -1.46695089e+01, -7.22850746e+01,
271 -1.65563055e+02, -3.10820425e+02, -4.70026655e+02,
272 -3.10820425e+02, -1.65563055e+02, -7.22850746e+01,
273 -1.46695089e+01, -3.37950835e+00, -6.91419168e-02,
274 3.57543979e+00, 1.89144704e+00],
275 [3.11841913e+00, 7.84024994e+00, 1.88495248e+00,
276 -7.69011009e+00, -2.71782400e+01, -1.04343326e+02,
277 -2.47561370e+02, -5.32959841e+02, -1.16529012e+03,
278 -5.32959841e+02, -2.47561370e+02, -1.04343326e+02,
279 -2.71782400e+01, -7.69011009e+00, 1.88495248e+00,
280 7.84024994e+00, 3.11841913e+00],
281 [2.74197956e+00, 4.73107997e+00, -9.48352966e-01,
282 -9.44822832e+00, -3.06477671e+01, -1.26788739e+02,
283 -3.22828411e+02, -8.47943472e+02, -3.87702420e+03,
284 -8.47943472e+02, -3.22828411e+02, -1.26788739e+02,
285 -3.06477671e+01, -9.44822832e+00, -9.48352966e-01,
286 4.73107997e+00, 2.74197956e+00],
287 [3.11841913e+00, 7.84024994e+00, 1.88495248e+00,
288 -7.69011009e+00, -2.71782400e+01, -1.04343326e+02,
289 -2.47561370e+02, -5.32959841e+02, -1.16529012e+03,
290 -5.32959841e+02, -2.47561370e+02, -1.04343326e+02,
291 -2.71782400e+01, -7.69011009e+00, 1.88495248e+00,
292 7.84024994e+00, 3.11841913e+00],
293 [1.89144704e+00, 3.57543979e+00, -6.91419168e-02,
294 -3.37950835e+00, -1.46695089e+01, -7.22850746e+01,
295 -1.65563055e+02, -3.10820425e+02, -4.70026655e+02,
296 -3.10820425e+02, -1.65563055e+02, -7.22850746e+01,
297 -1.46695089e+01, -3.37950835e+00, -6.91419168e-02,
298 3.57543979e+00, 1.89144704e+00],
299 [8.71929228e-01, 5.41025574e-01, 9.47560771e-01,
300 -5.75314708e-01, -7.46104027e+00, -4.42314481e+01,
301 -9.54126971e+01, -1.61603201e+02, -2.07520692e+02,
302 -1.61603201e+02, -9.54126971e+01, -4.42314481e+01,
303 -7.46104027e+00, -5.75314708e-01, 9.47560771e-01,
304 5.41025574e-01, 8.71929228e-01],
305 [1.05524430e+00, 1.71917897e+00, 1.73105590e+00,
306 -2.10088420e+00, -1.15118208e+01, -2.55007598e+01,
307 -4.73056159e+01, -6.97257685e+01, -8.09264433e+01,
308 -6.97257685e+01, -4.73056159e+01, -2.55007598e+01,
309 -1.15118208e+01, -2.10088420e+00, 1.73105590e+00,
310 1.71917897e+00, 1.05524430e+00],
311 [1.62986900e+00, 3.67851228e+00, 5.68645252e+00,
312 2.15342566e-01, -8.89937202e+00, -1.44739813e+01,
313 -2.98952660e+01, -4.37420817e+01, -4.83160958e+01,
314 -4.37420817e+01, -2.98952660e+01, -1.44739813e+01,
315 -8.89937202e+00, 2.15342566e-01, 5.68645252e+00,
316 3.67851228e+00, 1.62986900e+00],
317 [1.78571940e+00, 4.38918110e+00, 3.95098587e+00,
318 3.70961649e-01, -3.48151981e+00, -9.61567736e+00,
319 -1.78621172e+01, -2.32278872e+01, -2.31833727e+01,
320 -2.32278872e+01, -1.78621172e+01, -9.61567736e+00,
321 -3.48151981e+00, 3.70961649e-01, 3.95098587e+00,
322 4.38918110e+00, 1.78571940e+00],
323 [1.12382749e+00, 2.22609074e+00, 1.27877807e+00,
324 4.55434098e-01, -1.76842385e+00, -1.90046460e+00,
325 -8.10874526e+00, -1.20534899e+01, -1.48627948e+01,
326 -1.20534899e+01, -8.10874526e+00, -1.90046460e+00,
327 -1.76842385e+00, 4.55434098e-01, 1.27877807e+00,
328 2.22609074e+00, 1.12382749e+00],
329 [4.83499829e-01, 8.10171823e-01, 5.31096720e-01,
330 3.54369868e-02, -8.44782871e-01, -1.64614462+00,
331 -3.83933101e+00, -5.60243416e+00, -6.51691578e+00,
332 -5.60243416e+00, -3.83933101e+00, -1.64614462e+00,
333 -8.44782871e-01, 3.54369868e-02, 5.31096720e-01,
334 8.10171823e-01, 4.83499829e-01]]) * 1e-10
335 self.aN = np.array([[-9.672222587932733e-07,
336 -1.854684055704316e-07,
337 -6.907109379361769e-08,
338 -3.434510010399494e-08,
339 -1.9612893321972894e-08,
340 -1.2070595598419123e-08,
341 -7.76166581456265e-09,
342 -5.1290591302737845e-09],
343 [-1.6934304508907508e-07,
344 -1.1294085770257732e-07,
345 -5.556904377844212e-08,
346 -3.0404816410630946e-08,
347 -1.8115835808147195e-08,
348 -1.1400565570822897e-08,
349 -7.428030253054874e-09,
350 -4.950482905137934e-09],
351 [-2.642765928047855e-08,
352 -4.2973020270862075e-08,
353 -3.288284786525552e-08,
354 -2.2008013659494064e-08,
355 -1.4526432244373923e-08,
356 -9.685865982318977e-09,
357 -6.540042976345181e-09,
358 -4.463016792299867e-09],
359 [-7.890607506673887e-09,
360 -1.7268817838939093e-08,
361 -1.7571215335862644e-08,
362 -1.4190808979745564e-08,
363 -1.0525091813951477e-08,
364 -7.557427069020145e-09,
365 -5.3601997909453895e-09,
366 -3.785262753932024e-09],
367 [-3.19819383283163e-09,
368 -7.907188860704256e-09,
369 -9.431108045148581e-09,
370 -8.734433800366836e-09,
371 -7.177446016048238e-09,
372 -5.548891364163135e-09,
373 -4.151829718274712e-09,
374 -3.050079142680274e-09],
375 [-1.5264788341295506e-09,
376 -4.007588259092475e-09,
377 -5.237227746197923e-09,
378 -5.325369481025145e-09,
379 -4.74494467703471e-09,
380 -3.915820474757527e-09,
381 -3.084461007593805e-09,
382 -2.3593738079733816e-09],
383 [-8.028826757431073e-10,
384 -2.1802562891610477e-09,
385 -3.013496348411789e-09,
386 -3.2665536558956744e-09,
387 -3.0954503870565042e-09,
388 -2.6975377326297753e-09,
389 -2.225243828554418e-09,
390 -1.7689041096437489e-09],
391 [-4.494359608433485e-10,
392 -1.2463111600402407e-09,
393 -1.7870048664861338e-09,
394 -2.0261729415646362e-09,
395 -2.0119440116912244e-09,
396 -1.8328866609412017e-09,
397 -1.5737323964832772e-09,
398 -1.2957141987086006e-09]])
400 self.aE = np.array([[-6.777496868082258e-07,
401 -1.8380876162409285e-07,
402 -2.8018297221057767e-08,
403 -8.125248034789441e-09,
404 -3.2554136393266677e-09,
405 -1.5453475805766551e-09,
406 -8.104328349953093e-10,
407 -4.5288455522821904e-10],
408 [-1.6775243392445585e-07,
409 -1.0841387428867281e-07,
410 -4.332756930914473e-08,
411 -1.7516443134733466e-08,
412 -8.003213920226987e-09,
413 -4.046641342303334e-09,
414 -2.1977447699714067e-09,
415 -1.254855626337421e-09],
416 [-6.621450619914244e-08,
417 -5.399929099596603e-08,
418 -3.2620205287142746e-08,
419 -1.7612638233905413e-08,
420 -9.483506120439713e-09,
421 -5.269116790156857e-09,
422 -3.0311506399280317e-09,
423 -1.796806758939136e-09],
424 [-3.355107919271619e-08,
425 -2.983376538822224e-08,
426 -2.17864147917854e-08,
427 -1.4148255610378522e-08,
428 -8.742901288356405e-09,
429 -5.3404306089421096e-09,
430 -3.2783561039512352e-09,
431 -2.0341064582331245e-09],
432 [-1.93141015569346e-08,
433 -1.7872788139275774e-08,
434 -1.4392557125221263e-08,
435 -1.047378455931124e-08,
436 -7.165833862418497e-09,
437 -4.7471363350194744e-09,
438 -3.1007044557148345e-09,
439 -2.016814227704664e-09],
440 [-1.1934920937896483e-08,
441 -1.1282718164497807e-08,
442 -9.6076688011492e-09,
443 -7.516823026913443e-09,
444 -5.532490546520142e-09,
445 -3.911518551225112e-09,
446 -2.6981370339016726e-09,
447 -1.834970550556636e-09],
448 [-7.691827379621282e-09,
449 -7.364973173518246e-09,
450 -6.493299542676226e-09,
451 -5.331291165499783e-09,
452 -4.136757781512933e-09,
453 -3.077998169355289e-09,
454 -2.2233267681078023e-09,
455 -1.5738811344420306e-09],
456 [-5.089815731398664e-09,
457 -4.914168227172165e-09,
458 -4.434118887867189e-09,
459 -3.765280711029128e-09,
460 -3.0379496001714635e-09,
461 -2.3529138563027853e-09,
462 -1.76598680955912e-09,
463 -1.2947658524603495e-09]])
465 # Spline trap coefficients and the ctiCalibDict are all taken from a
466 # cti calibration measured from LSSTCam sensor R03_S12 during Run 5
467 # EO testing. These are the coefficients for the spline trap model
468 # used in the deferred charge calibration. The collection can be
469 # found in /repo/ir2: u/abrought/13144/cti (processed 3/4/2024).
470 self.splineTrapCoeffs = np.array([0.0, 28.1, 47.4, 56.4, 66.6, 78.6, 92.4, 109.4,
471 129.0, 151.9, 179.4, 211.9, 250.5, 296.2, 350.0,
472 413.5, 488.0, 576.0, 680.4, 753.0, 888.2, 1040.5,
473 1254.1, 1478.9, 1747.0, 2055.7, 2416.9, 2855.2,
474 3361.9, 3969.4, 4665.9, 5405.3, 6380.0, 7516.7,
475 8875.9, 10488.6, 12681.9, 14974.2, 17257.6, 20366.5,
476 24026.7, 28372.1, 33451.7, 39550.4, 46624.8, 55042.9,
477 64862.7, 76503.1, 90265.6, 106384.2, 0.0, 0.0, 0.0,
478 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.1,
479 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.1, 0.0,
480 0.1, 0.2, 0.6, 0.1, 0.0, 0.1, 0.0, 0.6, 0.3, 0.5, 0.8,
481 0.8, 1.5, 2.0, 1.8, 2.4, 2.6, 3.7, 5.0, 6.4, 8.4, 10.9,
482 14.5, 21.1, 28.9])
484 self.ctiCalibDict = {'inputGain': {'C:0,0': 1.5,
485 'C:0,1': 1.5,
486 'C:0,2': 1.5,
487 'C:0,3': 1.5,
488 'C:1,0': 1.5,
489 'C:1,1': 1.5,
490 'C:1,2': 1.5,
491 'C:1,3': 1.5},
492 'driftScale': {'C:0,0': 0.000127,
493 'C:0,1': 0.000137,
494 'C:0,2': 0.000138,
495 'C:0,3': 0.000147,
496 'C:1,0': 0.000147,
497 'C:1,1': 0.000122,
498 'C:1,2': 0.000123,
499 'C:1,3': 0.000116},
500 'decayTime': {'C:0,0': 2.30,
501 'C:0,1': 2.21,
502 'C:0,2': 2.28,
503 'C:0,3': 2.34,
504 'C:1,0': 2.30,
505 'C:1,1': 2.40,
506 'C:1,2': 2.51,
507 'C:1,3': 2.21},
508 'globalCti': {'C:0,0': 5.25e-07,
509 'C:0,1': 5.38e-07,
510 'C:0,2': 5.80e-07,
511 'C:0,3': 5.91e-07,
512 'C:1,0': 6.24e-07,
513 'C:1,1': 5.72e-07,
514 'C:1,2': 5.60e-07,
515 'C:1,3': 4.40e-07},
516 'signals': {'C:0,0': np.linspace(1.0e2, 1.0e5, 10),
517 'C:0,1': np.linspace(1.0e2, 1.0e5, 10),
518 'C:0,2': np.linspace(1.0e2, 1.0e5, 10),
519 'C:0,3': np.linspace(1.0e2, 1.0e5, 10),
520 'C:1,0': np.linspace(1.0e2, 1.0e5, 10),
521 'C:1,1': np.linspace(1.0e2, 1.0e5, 10),
522 'C:1,2': np.linspace(1.0e2, 1.0e5, 10),
523 'C:1,3': np.linspace(1.0e2, 1.0e5, 10)},
524 'serialEper': {'C:0,0': np.full(10, 5.25e-07),
525 'C:0,1': np.full(10, 5.38e-07),
526 'C:0,2': np.full(10, 5.80e-07),
527 'C:0,3': np.full(10, 5.91e-07),
528 'C:1,0': np.full(10, 6.24e-07),
529 'C:1,1': np.full(10, 5.72e-07),
530 'C:1,2': np.full(10, 5.60e-07),
531 'C:1,3': np.full(10, 4.40e-07)},
532 'parallelEper': {'C:0,0': np.full(10, 5.25e-07),
533 'C:0,1': np.full(10, 5.38e-07),
534 'C:0,2': np.full(10, 5.80e-07),
535 'C:0,3': np.full(10, 5.91e-07),
536 'C:1,0': np.full(10, 6.24e-07),
537 'C:1,1': np.full(10, 5.72e-07),
538 'C:1,2': np.full(10, 5.60e-07),
539 'C:1,3': np.full(10, 4.40e-07)},
540 'serialCtiTurnoff': {'C:0,0': 1.0e5,
541 'C:0,1': 1.0e5,
542 'C:0,2': 1.0e5,
543 'C:0,3': 1.0e5,
544 'C:1,0': 1.0e5,
545 'C:1,1': 1.0e5,
546 'C:1,2': 1.0e5,
547 'C:1,3': 1.0e5},
548 'parallelCtiTurnoff': {'C:0,0': 1.0e5,
549 'C:0,1': 1.0e5,
550 'C:0,2': 1.0e5,
551 'C:0,3': 1.0e5,
552 'C:1,0': 1.0e5,
553 'C:1,1': 1.0e5,
554 'C:1,2': 1.0e5,
555 'C:1,3': 1.0e5},
556 'serialCtiTurnoffSamplingErr': {'C:0,0': 1.0e3,
557 'C:0,1': 1.0e3,
558 'C:0,2': 1.0e3,
559 'C:0,3': 1.0e3,
560 'C:1,0': 1.0e3,
561 'C:1,1': 1.0e3,
562 'C:1,2': 1.0e3,
563 'C:1,3': 1.0e3},
564 'parallelCtiTurnoffSamplingErr': {'C:0,0': 1.0e3,
565 'C:0,1': 1.0e3,
566 'C:0,2': 1.0e3,
567 'C:0,3': 1.0e3,
568 'C:1,0': 1.0e3,
569 'C:1,1': 1.0e3,
570 'C:1,2': 1.0e3,
571 'C:1,3': 1.0e3},
572 'serialTraps': {'C:0,0': {'size': 20000.0,
573 'emissionTime': 0.4,
574 'pixel': 1,
575 'trap_type': 'spline',
576 'coeffs': self.splineTrapCoeffs},
577 'C:0,1': {'size': 20000.0,
578 'emissionTime': 0.4,
579 'pixel': 1,
580 'trap_type': 'spline',
581 'coeffs': self.splineTrapCoeffs},
582 'C:0,2': {'size': 20000.0,
583 'emissionTime': 0.4,
584 'pixel': 1,
585 'trap_type': 'spline',
586 'coeffs': self.splineTrapCoeffs},
587 'C:0,3': {'size': 20000.0,
588 'emissionTime': 0.4,
589 'pixel': 1,
590 'trap_type': 'spline',
591 'coeffs': self.splineTrapCoeffs},
592 'C:1,0': {'size': 20000.0,
593 'emissionTime': 0.4,
594 'pixel': 1,
595 'trap_type': 'spline',
596 'coeffs': self.splineTrapCoeffs},
597 'C:1,1': {'size': 20000.0,
598 'emissionTime': 0.4,
599 'pixel': 1,
600 'trap_type': 'spline',
601 'coeffs': self.splineTrapCoeffs},
602 'C:1,2': {'size': 20000.0,
603 'emissionTime': 0.4,
604 'pixel': 1,
605 'trap_type': 'spline',
606 'coeffs': self.splineTrapCoeffs},
607 'C:1,3': {'size': 20000.0,
608 'emissionTime': 0.4,
609 'pixel': 1,
610 'trap_type': 'spline',
611 'coeffs': self.splineTrapCoeffs}}}
613 self.deferredChargeCalib = self.makeDeferredChargeCalib()
615 # Cross-talk coeffs are defined in the parent class.
617 self.makeSubtask("assembleCcd")
619 def run(self):
620 """Generate a mock ISR product following LSSTCam ISR, and return it.
622 Returns
623 -------
624 image : `lsst.afw.image.Exposure`
625 Simulated ISR image with signals added.
626 dataProduct :
627 Simulated ISR data products.
628 None :
629 Returned if no valid configuration was found.
631 Raises
632 ------
633 RuntimeError
634 Raised if both doGenerateImage and doGenerateData are specified.
635 """
636 if self.config.doGenerateImage and self.config.doGenerateData:
637 raise RuntimeError("Only one of doGenerateImage and doGenerateData may be specified.")
638 elif self.config.doGenerateImage:
639 return self.makeImage()
640 elif self.config.doGenerateData:
641 return self.makeData()
642 else:
643 return None
645 def makeImage(self):
646 """Generate a simulated ISR LSST image.
648 Returns
649 -------
650 exposure : `lsst.afw.image.Exposure` or `dict`
651 Simulated ISR image data.
653 Notes
654 -----
655 This method constructs a "raw" data image.
656 """
657 exposure = self.getExposure()
659 # Set up random number generators for consistency of components,
660 # no matter the group that are configured.
661 rngSky = np.random.RandomState(seed=self.config.rngSeed + 1)
662 rngDark = np.random.RandomState(seed=self.config.rngSeed + 2)
663 rng2DBias = np.random.RandomState(seed=self.config.rngSeed + 3)
664 rngOverscan = np.random.RandomState(seed=self.config.rngSeed + 4)
665 rngReadNoise = np.random.RandomState(seed=self.config.rngSeed + 5)
666 rngBrighterFatter = galsim.BaseDeviate(self.config.rngSeed + 6)
668 # Create the linearizer if we will need it.
669 if self.config.doAddHighSignalNonlinearity:
670 linearizer = LinearizerMockLSST().run()
672 # We introduce effects as they happen from a source to the signal,
673 # so the effects go from electron to adu.
674 # The ISR steps will then correct these effects in the reverse order.
675 for idx, amp in enumerate(exposure.getDetector()):
676 # Get image bbox and data
677 bbox = None
678 if self.config.isTrimmed:
679 bbox = amp.getBBox()
680 bboxFull = bbox
681 else:
682 bbox = amp.getRawDataBBox()
683 bboxFull = amp.getRawBBox()
685 # This is the image data (excluding pre/overscans).
686 ampImageData = exposure.image[bbox]
687 # This is the full data (including pre/overscans if untrimmed).
688 ampFullData = exposure.image[bboxFull]
690 # Astrophysical signals are all in electron (e-).
691 # These are only applied to the imaging portion of the
692 # amplifier (ampImageData)
694 if self.config.doAddSky:
695 # The sky effects are in electron.
696 self.amplifierAddNoise(
697 ampImageData,
698 self.config.skyLevel,
699 np.sqrt(self.config.skyLevel),
700 rng=rngSky,
701 )
703 if self.config.doAddSource:
704 for sourceAmp, sourceFlux, sourceX, sourceY in zip(self.config.sourceAmp,
705 self.config.sourceFlux,
706 self.config.sourceX,
707 self.config.sourceY):
708 if idx == sourceAmp:
709 # The source flux is in electron.
710 self.amplifierAddSource(ampImageData, sourceFlux, sourceX, sourceY)
712 if self.config.doAddFringe:
713 # Fringes are added in electron.
714 self.amplifierAddFringe(amp,
715 ampImageData,
716 np.array(self.config.fringeScale),
717 x0=np.array(self.config.fringeX0),
718 y0=np.array(self.config.fringeY0))
720 if self.config.doAddFlat:
721 if self.config.calibMode:
722 # In case we are making a combined flat,
723 # add a non-zero signal so the mock flat can be multiplied
724 self.amplifierAddNoise(ampImageData, 1.0, 0.0)
725 # Multiply each amplifier by a Gaussian centered on u0 and v0
726 u0 = exposure.getDetector().getBBox().getDimensions().getX()/2.
727 v0 = exposure.getDetector().getBBox().getDimensions().getY()/2.
728 self.amplifierMultiplyFlat(amp, ampImageData, self.config.flatDrop, u0=u0, v0=v0)
730 # On-chip electronic effects.
732 # 1. Add bright defect(s).
733 if self.config.doAddBrightDefects:
734 defectList = self.makeDefectList(isTrimmed=self.config.isTrimmed)
736 for defect in defectList:
737 exposure.image[defect.getBBox()] = self.config.brightDefectLevel
739 for idx, amp in enumerate(exposure.getDetector()):
740 # Get image bbox and data
741 bbox = None
742 if self.config.isTrimmed:
743 bbox = amp.getBBox()
744 bboxFull = bbox
745 else:
746 bbox = amp.getRawDataBBox()
747 bboxFull = amp.getRawBBox()
749 # This is the image data (excluding pre/overscans).
750 ampImageData = exposure.image[bbox]
751 # This is the full data (including pre/overscans if untrimmed).
752 ampFullData = exposure.image[bboxFull]
754 # 2. Add dark current (electron) to imaging portion of the amp.
755 if self.config.doAddDark or self.config.doAddDarkNoiseOnly:
756 if self.config.doAddDarkNoiseOnly:
757 darkLevel = 0.0
758 else:
759 darkLevel = self.config.darkRate * self.config.darkTime
760 if self.config.calibMode:
761 darkNoise = 0.0
762 else:
763 darkNoise = np.sqrt(self.config.darkRate * self.config.darkTime)
765 self.amplifierAddNoise(ampImageData, darkLevel, darkNoise, rng=rngDark)
767 # 3. Add BF effect (electron) to imaging portion of the amp.
768 if self.config.doAddBrighterFatter is True:
769 self.amplifierAddBrighterFatter(
770 ampImageData,
771 rngBrighterFatter,
772 self.config.bfStrength,
773 self.config.nRecalc,
774 )
776 # 4. Add serial CTI (electron) to amplifier (imaging + overscan).
777 if self.config.doAddDeferredCharge:
778 # Get the free charge area for the amplifier.
779 self.amplifierAddDeferredCharge(exposure, amp)
781 # 5. Add 2D bias residual (electron) to imaging portion of the amp.
782 if self.config.doAdd2DBias:
783 # For now we use an unstructured noise field to add some
784 # consistent 2D bias residual that can be subtracted. In
785 # the future this can be made into a warm corner (for example).
786 self.amplifierAddNoise(
787 ampImageData,
788 0.0,
789 self.config.noise2DBias,
790 rng=rng2DBias,
791 )
793 # 6. Add clock-injected offset (electron) to amplifer
794 # (imaging + overscan).
795 # This is just an offset that will be crosstalked and modified by
796 # the gain, and does not have a noise associated with it.
797 if self.config.doAddClockInjectedOffset:
798 self.amplifierAddNoise(
799 ampFullData,
800 self.config.clockInjectedOffsetLevel,
801 0.0,
802 )
804 # 7./8. Add serial and parallel overscan slopes (electron)
805 # (imaging + overscan)
806 if (self.config.doAddParallelOverscanRamp or self.config.doAddSerialOverscanRamp) and \
807 not self.config.isTrimmed:
809 if self.config.doAddParallelOverscanRamp:
810 # Apply gradient along the X axis.
811 self.amplifierAddXGradient(ampFullData, -1.0 * self.config.overscanScale,
812 1.0 * self.config.overscanScale)
814 if self.config.doAddSerialOverscanRamp:
815 # Apply the gradient along the Y axis.
816 self.amplifierAddYGradient(ampFullData, -1.0 * self.config.overscanScale,
817 1.0 * self.config.overscanScale)
819 # 9. Add non-linearity (electron) to amplifier
820 # (imaging + overscan).
821 if self.config.doAddHighSignalNonlinearity:
822 # The linearizer coefficients come from makeLinearizer().
823 if linearizer.linearityType[amp.getName()] != "Spline":
824 raise RuntimeError("IsrMockLSST only supports spline non-linearity.")
826 coeffs = linearizer.linearityCoeffs[amp.getName()]
827 centers, values = np.split(coeffs, 2)
829 # This is an application of high signal non-linearity, so we
830 # set the lower values to 0.0 (this cut is arbitrary).
831 values[centers < self.config.highSignalNonlinearityThreshold] = 0.0
833 # The linearizer is units of adu, so convert to electron
834 centers *= self.config.gainDict[amp.getName()]
835 values *= self.config.gainDict[amp.getName()]
837 # Note that the linearity spline is in "overscan subtracted"
838 # units so needs to be applied without the clock-injected
839 # offset.
840 self.amplifierAddNonlinearity(
841 ampFullData,
842 centers,
843 values,
844 self.config.clockInjectedOffsetLevel if self.config.doAddClockInjectedOffset else 0.0,
845 )
847 # 10. Add read noise (electron) to the amplifier
848 # (imaging + overscan).
849 # Unsure if this should be before or after crosstalk.
850 # Probably some of both; hopefully doesn't matter.
851 if not self.config.calibMode:
852 # Add read noise to the imaging region.
853 self.amplifierAddNoise(
854 ampImageData,
855 0.0,
856 self.config.readNoise,
857 rng=rngReadNoise,
858 )
860 # If not trimmed, add to the overscan regions.
861 if not self.config.isTrimmed:
862 parallelOverscanBBox = amp.getRawParallelOverscanBBox()
863 parallelOverscanData = exposure.image[parallelOverscanBBox]
865 serialOverscanBBox = self.getFullSerialOverscanBBox(amp)
866 serialOverscanData = exposure.image[serialOverscanBBox]
868 # Add read noise of mean 0
869 # to the parallel and serial overscan regions.
870 self.amplifierAddNoise(
871 parallelOverscanData,
872 0.0,
873 self.config.readNoise,
874 rng=rngOverscan,
875 )
876 self.amplifierAddNoise(
877 serialOverscanData,
878 0.0,
879 self.config.readNoise,
880 rng=rngOverscan,
881 )
883 # 7b. Add bad column to the parallel overscan region.
884 if self.config.doAddBadParallelOverscanColumn and not self.config.isTrimmed:
885 # We want to place this right above the defect, to simulate
886 # bleeding into the parallel overscan region.
887 amp = exposure.getDetector()[2]
888 parBBox = amp.getRawParallelOverscanBBox()
889 bboxBad = geom.Box2I(
890 corner=geom.Point2I(50, parBBox.getMinY()),
891 dimensions=geom.Extent2I(1, parBBox.getHeight()),
892 )
893 exposure[bboxBad].image.array[:, :] += self.config.badParallelOverscanColumnLevel
895 if self.config.doAddBadParallelOverscanColumnNeighbors:
896 for neighbor in [49, 51]:
897 bboxBad = geom.Box2I(
898 corner=geom.Point2I(neighbor, parBBox.getMinY()),
899 dimensions=geom.Extent2I(1, parBBox.getHeight()),
900 )
901 exposure[bboxBad].image.array[:, :] += self.config.badParallelOverscanColumnNeighborsLevel
903 # 11. Add crosstalk (electron) to all the amplifiers
904 # (imaging + overscan).
905 if self.config.doAddCrosstalk:
906 ctCalib = CrosstalkCalib()
907 # We use the regular subtractCrosstalk code but with a negative
908 # sign on the crosstalk coefficients so it adds instead of
909 # subtracts. We only apply the signal plane (ignoreVariance,
910 # subtrahendMasking) with a very large pixel to mask to ensure
911 # no crosstalk mask bits are set.
912 ctCalib.subtractCrosstalk(
913 exposure,
914 crosstalkCoeffs=-1*self.crosstalkCoeffs,
915 doSubtrahendMasking=True,
916 minPixelToMask=np.inf,
917 ignoreVariance=True,
918 fullAmplifier=True,
919 )
921 for amp in exposure.getDetector():
922 # Get image bbox and data (again).
923 bbox = None
924 if self.config.isTrimmed:
925 bbox = amp.getBBox()
926 bboxFull = bbox
927 else:
928 bbox = amp.getRawDataBBox()
929 bboxFull = amp.getRawBBox()
931 # This is the image data (excluding pre/overscans).
932 ampImageData = exposure.image[bbox]
933 # This is the full data (including pre/overscans if untrimmed).
934 ampFullData = exposure.image[bboxFull]
936 # 12. Gain un-normalize (from electron to floating point adu)
937 if self.config.doApplyGain:
938 gain = self.config.gainDict.get(amp.getName(), self.config.gain)
939 self.applyGain(ampFullData, gain)
941 # 13. Add overall bias level (adu) to the amplifier
942 # (imaging + overscan)
943 if self.config.doAddBias:
944 self.addBiasLevel(ampFullData, self.config.biasLevel)
946 # 14. Round/Truncate to integers (adu)
947 if self.config.doRoundAdu:
948 self.roundADU(ampFullData)
950 # Add units metadata to calibrations.
951 if self.config.calibMode:
952 if self.config.doApplyGain:
953 exposure.metadata["LSST ISR UNITS"] = "adu"
954 else:
955 exposure.metadata["LSST ISR UNITS"] = "electron"
957 # Add a variance plane appropriate for a calibration frame.
958 # We take the absolute value for biases which have no signal.
959 exposure.variance.array[:, :] = np.abs(np.median(exposure.image.array)/10.)
961 exposure.metadata["BSSVBS"] = 50.0
962 exposure.metadata["HVBIAS"] = "ON"
964 if self.config.doGenerateAmpDict:
965 expDict = dict()
966 for amp in exposure.getDetector():
967 expDict[amp.getName()] = exposure
968 return expDict
969 else:
970 return exposure
972 def addBiasLevel(self, ampData, biasLevel):
973 """Add bias level to an amplifier's image data.
975 Parameters
976 ----------
977 ampData : `lsst.afw.image.ImageF`
978 Amplifier image to operate on.
979 biasLevel : `float`
980 Bias level to be added to the image.
981 """
982 ampArr = ampData.array
983 ampArr[:] = ampArr[:] + biasLevel
985 def makeDefectList(self, isTrimmed=True):
986 """Generate a simple defect list.
988 Parameters
989 ----------
990 isTrimmed : `bool`, optional
991 Return defects in trimmed coordinates?
993 Returns
994 -------
995 defectList : `lsst.meas.algorithms.Defects`
996 Simulated defect list
997 """
998 defectBoxesUntrimmed = [
999 geom.Box2I(
1000 geom.Point2I(50, 118),
1001 geom.Extent2I(1, 51),
1002 ),
1003 ]
1005 if not isTrimmed:
1006 return Defects(defectBoxesUntrimmed)
1008 # If trimmed, we need to convert.
1009 tempExp = self.getExposure(isTrimmed=False)
1010 tempExp.image.array[:, :] = 0.0
1011 for bbox in defectBoxesUntrimmed:
1012 tempExp.image[bbox] = 1.0
1014 assembledExp = self.assembleCcd.assembleCcd(tempExp)
1016 # Use thresholding code to find defect footprints/boxes.
1017 threshold = afwDetection.createThreshold(1.0, "value", polarity=True)
1018 footprintSet = afwDetection.FootprintSet(assembledExp.image, threshold)
1020 return Defects.fromFootprintList(footprintSet.getFootprints())
1022 def makeBfKernel(self):
1023 """Generate a simple simulated brighter-fatter kernel.
1024 Returns
1025 -------
1026 kernel : `lsst.ip.isr.BrighterFatterKernel`
1027 Simulated brighter-fatter kernel.
1028 """
1029 bfkArray = super().makeBfKernel()
1030 bfKernelObject = BrighterFatterKernel()
1031 bfKernelObject.level = 'AMP'
1032 bfKernelObject.gain = self.config.gainDict
1034 for amp in self.getExposure().getDetector():
1035 # Kernel must be in (y,x) orientation
1036 bfKernelObject.ampKernels[amp.getName()] = bfkArray.T
1037 bfKernelObject.valid[amp.getName()] = True
1039 return bfKernelObject
1041 def makeElectrostaticBf(self):
1042 """Generate a simple simulated electrostatic
1043 brighter-fatter calibration.
1045 Returns
1046 -------
1047 kernel : `lsst.ip.isr.ElectrostaticBrighterFatterDistortionMatrix`
1048 Simulated brighter-fatter kernel.
1049 """
1050 aN, aS, aE, aW = super().makeElectrostaticBf()
1051 electroBfDistortionMatrix = ElectrostaticBrighterFatterDistortionMatrix(
1052 inputRange=self.aN.shape[0],
1053 fitRange=self.aN.shape[0],
1054 )
1055 electroBfDistortionMatrix.aN = self.aN
1056 electroBfDistortionMatrix.aS = self.aN
1057 electroBfDistortionMatrix.aE = self.aE
1058 electroBfDistortionMatrix.aW = self.aE # Assume it is symmetric
1059 electroBfDistortionMatrix.gain = self.config.gainDict
1061 return electroBfDistortionMatrix
1063 def makeDeferredChargeCalib(self):
1064 """Generate a CTI calibration.
1066 Returns
1067 -------
1068 cti : `lsst.ip.isr.deferredCharge.DeferredChargeCalib`
1069 Simulated deferred charge calibration.
1070 """
1072 metadataDict = {'metadata': PropertyList()}
1073 metadataDict['metadata'].add(name="OBSTYPE", value="CTI")
1074 metadataDict['metadata'].add(name="CALIBCLS",
1075 value="lsst.ip.isr.deferredCharge.DeferredChargeCalib")
1076 self.ctiCalibDict = {**metadataDict, **self.ctiCalibDict}
1077 deferredChargeCalib = DeferredChargeCalib()
1078 self.cti = deferredChargeCalib.fromDict(self.ctiCalibDict)
1080 return self.cti
1082 def amplifierAddBrighterFatter(self, ampImageData, rng, bfStrength, nRecalc):
1083 """Add brighter fatter effect and/or diffusion to the image.
1084 Parameters
1085 ----------
1086 ampImageData : `lsst.afw.image.ImageF`
1087 Amplifier image to operate on.
1088 rng : `galsim.BaseDeviate`
1089 Random number generator.
1090 bfStrength : `float`
1091 Scaling parameter of the brighter fatter effect (nominally = 1)
1092 nRecalc: 'int'
1093 The number of electrons to accumulate before recalculating the
1094 distortion of the pixel shapes.
1095 """
1097 incidentImage = galsim.Image(ampImageData.array, scale=1)
1098 measuredImage = galsim.ImageF(
1099 ampImageData.array.shape[1],
1100 ampImageData.array.shape[0],
1101 scale=1,
1102 )
1103 photons = galsim.PhotonArray.makeFromImage(incidentImage)
1105 sensorModel = galsim.SiliconSensor(
1106 strength=bfStrength,
1107 rng=rng,
1108 diffusion_factor=0.0,
1109 nrecalc=nRecalc,
1110 )
1112 totalFluxAdded = sensorModel.accumulate(photons, measuredImage)
1113 ampImageData.array = measuredImage.array
1115 return totalFluxAdded
1117 def amplifierAddDeferredCharge(self, exposure, amp):
1118 """Add serial CTI to the amplifier data.
1120 Parameters
1121 ----------
1122 exposure : `lsst.afw.image.ExposureF`
1123 The exposure object containing the amplifier
1124 to apply deferred charge to.
1125 amp : `lsst.afw.image.Amplifier`
1126 The amplifier object (contains geometry info).
1127 """
1128 # Get the amplifier's geometry parameters.
1129 # When adding deferred charge, we have already assured that
1130 # isTrimmed is False. Therefore we want to make sure that we
1131 # get the RawDataBBox.
1132 readoutCorner = amp.getReadoutCorner()
1133 prescanWidth = amp.getRawHorizontalPrescanBBox().getWidth()
1134 serialOverscanWidth = amp.getRawHorizontalOverscanBBox().getWidth()
1135 parallelOverscanWidth = amp.getRawVerticalOverscanBBox().getHeight()
1136 bboxFreeCharge = amp.getRawDataBBox()
1137 bboxFreeCharge = bboxFreeCharge.expandedTo(amp.getRawHorizontalOverscanBBox())
1138 bboxFreeCharge = bboxFreeCharge.expandedTo(amp.getRawHorizontalPrescanBBox())
1139 bboxFreeCharge = bboxFreeCharge.expandedTo(amp.getRawVerticalOverscanBBox())
1141 ampFreeChargeData = exposure.image[bboxFreeCharge]
1142 ampImageData = exposure.image[amp.getRawDataBBox()]
1144 # Get the deferred charge parameters for this amplifier.
1145 cti = self.deferredChargeCalib.globalCti[amp.getName()]
1146 traps = self.deferredChargeCalib.serialTraps[amp.getName()]
1147 driftScale = self.deferredChargeCalib.driftScale[amp.getName()]
1148 decayTime = self.deferredChargeCalib.decayTime[amp.getName()]
1150 # Create a fake amplifier object that contains some deferred charge
1151 # paramters.
1152 floatingOutputAmplifier = FloatingOutputAmplifier(
1153 gain=1.0, # Everything is already in electrons.
1154 scale=driftScale,
1155 decay_time=decayTime,
1156 noise=0.0,
1157 offset=0.0,
1158 )
1160 def flipImage(arr, readoutCorner):
1161 # Flip an array so that the readout corner is in
1162 # the lower left.
1163 if readoutCorner == ReadoutCorner.LR:
1164 return np.fliplr(arr)
1165 elif readoutCorner == ReadoutCorner.UR:
1166 return np.fliplr(np.flipud(arr))
1167 elif readoutCorner == ReadoutCorner.UL:
1168 return np.flipud(arr)
1169 else:
1170 pass
1172 return arr
1174 # The algorithm expects that the readout corner is in
1175 # the lower left corner. Flip it to be so:
1176 imageData = ampImageData.array
1177 imageData = flipImage(imageData, readoutCorner)
1179 # Simulate the amplifier.
1180 ampSim = SegmentSimulator(
1181 imarr=imageData,
1182 prescan_width=prescanWidth,
1183 output_amplifier=floatingOutputAmplifier,
1184 cti=cti,
1185 traps=traps,
1186 )
1188 # Simulate deferred charge!
1189 # Note that the readout() method uses the image region data and the
1190 # overscan dimensions provided as input parameters. It then creates
1191 # overscan nd adds it to the image data to create the raw image.
1192 result = ampSim.readout(
1193 serial_overscan_width=serialOverscanWidth,
1194 parallel_overscan_width=parallelOverscanWidth,
1195 )
1197 # Flip the image back to the original orientation.
1198 result = flipImage(result, readoutCorner)
1200 # Set the image with the deferred charge added.
1201 ampFreeChargeData.array[:, :] = result
1203 def makeLinearizer(self):
1204 # docstring inherited.
1206 # The linearizer has units of adu.
1207 nNodes = 10
1208 # Set this to just above the mock saturation (adu)
1209 maxADU = 101_000
1210 nonLinSplineNodes = np.linspace(0, maxADU, nNodes)
1211 # These values come from cp_pipe/tests/test_linearity.py and
1212 # are based on a test fit to LSSTCam data, run 7193D, detector 22,
1213 # amp C00.
1214 nonLinSplineValues = np.array(
1215 [0.0, -8.87, 1.46, 1.69, -6.92, -68.23, -78.01, -11.56, 80.26, 185.01]
1216 )
1218 if self.config.doAddHighSignalNonlinearity and not self.config.doAddLowSignalNonlinearity:
1219 nonLinSplineValues[nonLinSplineNodes < self.config.highSignalNonlinearityThreshold] = 0.0
1220 elif self.config.doAddLowSignalNonlinearity:
1221 raise NotImplementedError("Low signal non-linearity is not implemented.")
1223 exp = self.getExposure()
1224 detector = exp.getDetector()
1226 linearizer = Linearizer(detector=detector)
1227 linearizer.updateMetadataFromExposures([exp])
1229 # We need to set override by hand because we are constructing a
1230 # linearizer manually and not from a serialized object.
1231 linearizer.override = True
1232 linearizer.hasLinearity = True
1233 linearizer.validate()
1234 linearizer.updateMetadata(camera=self.getCamera(), detector=detector, filterName='NONE')
1235 linearizer.updateMetadata(setDate=True, setCalibId=True)
1237 for amp in detector:
1238 ampName = amp.getName()
1239 linearizer.linearityType[ampName] = "Spline"
1240 linearizer.linearityCoeffs[ampName] = np.concatenate([nonLinSplineNodes, nonLinSplineValues])
1241 # We need to specify the raw bbox here.
1242 linearizer.linearityBBox[ampName] = amp.getRawBBox()
1244 return linearizer
1246 def amplifierAddNonlinearity(self, ampData, centers, values, offset):
1247 """Add non-linearity to amplifier data.
1249 Parameters
1250 ----------
1251 ampData : `lsst.afw.image.ImageF`
1252 Amplifier image to operate on.
1253 centers : `np.ndarray`
1254 Spline nodes.
1255 values : `np.ndarray`
1256 Spline values.
1257 offset : `float`
1258 Offset zero-point between linearizer (internal vs external).
1259 """
1260 # I'm not sure what to do about negative values...
1262 # Note that we are using the afw AKIMA_SPLINE to offset the
1263 # data but using the equivalent but faster scipy Akima1DInterpolator to
1264 # correct the data.
1265 spl = afwMath.makeInterpolate(
1266 centers,
1267 values,
1268 afwMath.stringToInterpStyle("AKIMA_SPLINE"),
1269 )
1271 delta = np.asarray(spl.interpolate(ampData.array.ravel() - offset))
1273 ampData.array[:, :] += delta.reshape(ampData.array.shape)
1275 def amplifierMultiplyFlat(self, amp, ampData, fracDrop, u0=100.0, v0=100.0):
1276 """Multiply an amplifier's image data by a flat-like pattern.
1278 Parameters
1279 ----------
1280 amp : `lsst.afw.ampInfo.AmpInfoRecord`
1281 Amplifier to operate on. Needed for amp<->exp coordinate
1282 transforms.
1283 ampData : `lsst.afw.image.ImageF`
1284 Amplifier image to operate on.
1285 fracDrop : `float`
1286 Fractional drop from center to edge of detector along x-axis.
1287 u0 : `float`
1288 Peak location in detector coordinates.
1289 v0 : `float`
1290 Peak location in detector coordinates.
1291 """
1292 if fracDrop >= 1.0:
1293 raise RuntimeError("Flat fractional drop cannot be greater than 1.0")
1295 sigma = u0 / np.sqrt(2.0 * fracDrop)
1297 for x in range(0, ampData.getDimensions().getX()):
1298 for y in range(0, ampData.getDimensions().getY()):
1299 (u, v) = self.localCoordToExpCoord(amp, x, y)
1300 f = np.exp(-0.5 * ((u - u0)**2 + (v - v0)**2) / sigma**2)
1301 ampData.array[y][x] = (ampData.array[y][x] * f)
1303 def applyGain(self, ampData, gain):
1304 """Apply gain to the amplifier's data.
1305 This method divides the data by the gain
1306 because the mocks need to convert the data in electron to adu,
1307 so it does the inverse operation to applyGains in isrFunctions.
1309 Parameters
1310 ----------
1311 ampData : `lsst.afw.image.ImageF`
1312 Amplifier image to operate on.
1313 gain : `float`
1314 Gain value in electron/adu.
1315 """
1316 ampArr = ampData.array
1317 ampArr[:] = ampArr[:] / gain
1319 def roundADU(self, ampData):
1320 """Round adu to nearest integer.
1322 Parameters
1323 ----------
1324 ampData : `lsst.afw.image.ImageF`
1325 Amplifier image to operate on.
1326 """
1327 ampArr = ampData.array
1328 ampArr[:] = np.around(ampArr)
1330 def amplifierAddXGradient(self, ampData, start, end):
1331 """Add a x-axis linear gradient to an amplifier's image data.
1333 This method operates in the amplifier coordinate frame.
1335 Parameters
1336 ----------
1337 ampData : `lsst.afw.image.ImageF`
1338 Amplifier image to operate on.
1339 start : `float`
1340 Start value of the gradient (at x=0).
1341 end : `float`
1342 End value of the gradient (at x=xmax).
1343 """
1344 nPixX = ampData.getDimensions().getX()
1345 ampArr = ampData.array
1346 ampArr[:] = ampArr[:] + (np.interp(range(nPixX), (0, nPixX - 1), (start, end)).reshape(1, nPixX)
1347 + np.zeros(ampData.getDimensions()).transpose())
1349 def getFullSerialOverscanBBox(self, amp):
1350 """Get the full serial overscan bounding box from an amplifier.
1352 This includes the serial/parallel overscan region.
1354 Parameters
1355 ----------
1356 amp : `lsst.afw.ampInfo.AmpInfoRecord`
1357 Amplifier to operate on.
1359 Returns
1360 -------
1361 bbox : `lsst.geom.Box2I`
1362 """
1363 # This only works for untrimmed data.
1364 bbox = amp.getRawDataBBox()
1366 parallelOverscanBBox = amp.getRawParallelOverscanBBox()
1367 grownImageBBox = bbox.expandedTo(parallelOverscanBBox)
1369 serialOverscanBBox = amp.getRawSerialOverscanBBox()
1370 # Extend the serial overscan bbox to include corners
1371 serialOverscanBBox = geom.Box2I(
1372 geom.Point2I(serialOverscanBBox.getMinX(),
1373 grownImageBBox.getMinY()),
1374 geom.Extent2I(serialOverscanBBox.getWidth(),
1375 grownImageBBox.getHeight()))
1377 return serialOverscanBBox
1380class RawMockLSST(IsrMockLSST):
1381 """Generate a raw exposure suitable for ISR.
1382 """
1383 def __init__(self, **kwargs):
1384 super().__init__(**kwargs)
1385 self.config.isTrimmed = False
1386 self.config.doGenerateImage = True
1387 self.config.doGenerateAmpDict = False
1389 # Add astro effects
1390 self.config.doAddSky = True
1391 self.config.doAddSource = True
1393 # Add optical effects
1394 self.config.doAddFringe = True
1396 # Add instrument effects
1397 self.config.doAddParallelOverscanRamp = True
1398 self.config.doAddSerialOverscanRamp = True
1399 self.config.doAddCrosstalk = True
1400 self.config.doAddBias = True
1401 self.config.doAddDark = True
1403 self.config.doAddFlat = True
1406class TrimmedRawMockLSST(RawMockLSST):
1407 """Generate a trimmed raw exposure.
1408 """
1409 def __init__(self, **kwargs):
1410 super().__init__(**kwargs)
1411 self.config.isTrimmed = True
1412 self.config.doAddParallelOverscanRamp = False
1413 self.config.doAddSerialOverscanRamp = False
1416class CalibratedRawMockLSST(RawMockLSST):
1417 """Generate a trimmed raw exposure.
1419 This represents a "truth" image that can be compared to a
1420 post-ISR cleaned image.
1421 """
1422 def __init__(self, **kwargs):
1423 super().__init__(**kwargs)
1424 self.config.isTrimmed = True
1425 self.config.doGenerateImage = True
1427 self.config.doAddSky = True
1428 self.config.doAddSource = True
1430 self.config.doAddFringe = True
1432 self.config.doAddParallelOverscanRamp = False
1433 self.config.doAddSerialOverscanRamp = False
1434 self.config.doAddCrosstalk = False
1435 self.config.doAddBias = False
1436 self.config.doAdd2DBias = False
1437 self.config.doAddDark = False
1438 self.config.doApplyGain = False
1439 self.config.doAddFlat = False
1440 self.config.doAddClockInjectedOffset = False
1442 self.config.biasLevel = 0.0
1443 # Assume combined calibrations are made with 16 inputs.
1444 self.config.readNoise *= 0.25
1446 self.config.doRoundAdu = False
1449class ReferenceMockLSST(IsrMockLSST):
1450 """Parent class for those that make reference calibrations.
1451 """
1452 def __init__(self, **kwargs):
1453 # If we want the calibration in adu units, we need to apply
1454 # the gain. Default is electron units, so do not apply the gain.
1455 doApplyGain = kwargs.pop("adu", False)
1457 super().__init__(**kwargs)
1458 self.config.isTrimmed = True
1459 self.config.doGenerateImage = True
1461 self.config.calibMode = True
1463 self.config.doAddSky = False
1464 self.config.doAddSource = False
1466 self.config.doAddFringe = False
1468 self.config.doAddParallelOverscanRamp = False
1469 self.config.doAddSerialOverscanRamp = False
1470 self.config.doAddCrosstalk = False
1471 self.config.doAddBias = False
1472 self.config.doAdd2DBias = False
1473 self.config.doAddDark = False
1474 self.config.doApplyGain = doApplyGain
1475 self.config.doAddFlat = False
1476 self.config.doAddClockInjectedOffset = False
1478 # Reference calibrations are not integerized.
1479 self.config.doRoundAdu = False
1482# Classes to generate calibration products mocks.
1483class DarkMockLSST(ReferenceMockLSST):
1484 """Simulated reference dark calibration.
1485 """
1486 def __init__(self, **kwargs):
1487 super().__init__(**kwargs)
1488 self.config.doAddDark = True
1489 self.config.darkTime = 1.0
1492class BiasMockLSST(ReferenceMockLSST):
1493 """Simulated combined bias calibration.
1494 """
1495 def __init__(self, **kwargs):
1496 super().__init__(**kwargs)
1497 # This is a "2D bias residual" frame which has only
1498 # the 2D bias in it.
1499 self.config.doAdd2DBias = True
1502class FlatMockLSST(ReferenceMockLSST):
1503 """Simulated reference flat calibration.
1504 """
1505 def __init__(self, **kwargs):
1506 super().__init__(**kwargs)
1507 self.config.doAddFlat = True
1510class FringeMockLSST(ReferenceMockLSST):
1511 """Simulated reference fringe calibration.
1512 """
1513 def __init__(self, **kwargs):
1514 super().__init__(**kwargs)
1515 self.config.doAddFringe = True
1518class BfKernelMockLSST(IsrMockLSST):
1519 """Simulated brighter-fatter kernel.
1520 """
1521 def __init__(self, **kwargs):
1522 super().__init__(**kwargs)
1523 self.config.doGenerateImage = False
1524 self.config.doGenerateData = True
1526 self.config.doBrighterFatter = True
1527 self.config.brighterFatterCalibType = "KERNEL"
1528 self.config.doDefects = False
1529 self.config.doCrosstalkCoeffs = False
1530 self.config.doTransmissionCurve = False
1531 self.config.doLinearizer = False
1534class ElectrostaticBfMockLSST(IsrMockLSST):
1535 """Simulated electrostatic brighter-fatter
1536 calibration.
1537 """
1538 def __init__(self, **kwargs):
1539 super().__init__(**kwargs)
1540 self.config.doGenerateImage = False
1541 self.config.doGenerateData = True
1543 self.config.doBrighterFatter = True
1544 self.config.brighterFatterCalibType = "ELECTROSTATIC"
1545 self.config.doDefects = False
1546 self.config.doCrosstalkCoeffs = False
1547 self.config.doTransmissionCurve = False
1548 self.config.doLinearizer = False
1551class DeferredChargeMockLSST(IsrMockLSST):
1552 """Simulated deferred charge calibration.
1553 """
1554 def __init__(self, **kwargs):
1555 super().__init__(**kwargs)
1556 self.config.doGenerateImage = False
1557 self.config.doGenerateData = True
1558 self.config.doDeferredCharge = True
1559 self.config.doDefects = False
1560 self.config.doCrosstalkCoeffs = False
1561 self.config.doTransmissionCurve = False
1564class DefectMockLSST(IsrMockLSST):
1565 """Simulated defect list.
1566 """
1567 def __init__(self, **kwargs):
1568 super().__init__(**kwargs)
1569 self.config.doGenerateImage = False
1570 self.config.doGenerateData = True
1572 self.config.doBrighterFatter = False
1573 self.config.doDefects = True
1574 self.config.doCrosstalkCoeffs = False
1575 self.config.doTransmissionCurve = False
1576 self.config.doLinearizer = False
1579class CrosstalkCoeffMockLSST(IsrMockLSST):
1580 """Simulated crosstalk coefficient matrix.
1581 """
1582 def __init__(self, **kwargs):
1583 super().__init__(**kwargs)
1584 self.config.doGenerateImage = False
1585 self.config.doGenerateData = True
1587 self.config.doBrighterFatter = False
1588 self.config.doDefects = False
1589 self.config.doCrosstalkCoeffs = True
1590 self.config.doTransmissionCurve = False
1591 self.config.doLinearizer = False
1594class LinearizerMockLSST(IsrMockLSST):
1595 """Simulated linearizer.
1596 """
1597 def __init__(self, **kwargs):
1598 super().__init__(**kwargs)
1599 self.config.doGenerateImage = False
1600 self.config.doGenerateData = True
1602 self.config.doBrighterFatter = False
1603 self.config.doDefects = False
1604 self.config.doCrosstalkCoeffs = False
1605 self.config.doTransmissionCurve = False
1606 self.config.doLinearizer = True
1609class TransmissionMockLSST(IsrMockLSST):
1610 """Simulated transmission curve.
1611 """
1612 def __init__(self, **kwargs):
1613 super().__init__(**kwargs)
1614 self.config.doGenerateImage = False
1615 self.config.doGenerateData = True
1617 self.config.doBrighterFatter = False
1618 self.config.doDefects = False
1619 self.config.doCrosstalkCoeffs = False
1620 self.config.doTransmissionCurve = True
1621 self.config.doLinearizer = False