Coverage for python/lsst/ip/isr/isrMockLSST.py: 16%

465 statements  

« prev     ^ index     » next       coverage.py v7.11.0, created at 2025-11-04 02:08 -0800

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__ = ["IsrMockLSSTConfig", "IsrMockLSST", "RawMockLSST", 

23 "CalibratedRawMockLSST", "ReferenceMockLSST", 

24 "BiasMockLSST", "DarkMockLSST", "FlatMockLSST", "FringeMockLSST", 

25 "BfKernelMockLSST", "DefectMockLSST", "CrosstalkCoeffMockLSST", 

26 "TransmissionMockLSST"] 

27 

28import numpy as np 

29import galsim 

30 

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 .deferredCharge import ( 

44 SegmentSimulator, 

45 FloatingOutputAmplifier, 

46 DeferredChargeCalib 

47) 

48 

49 

50class IsrMockLSSTConfig(IsrMockConfig): 

51 """Configuration parameters for isrMockLSST. 

52 """ 

53 # Detector parameters and "Exposure" parameters, 

54 # mostly inherited from IsrMockConfig. 

55 isLsstLike = pexConfig.Field( 

56 dtype=bool, 

57 default=True, 

58 doc="If True, products have one raw image per amplifier, otherwise, one raw image per detector.", 

59 ) 

60 calibMode = pexConfig.Field( 

61 dtype=bool, 

62 default=False, 

63 doc="Set to true to produce mock calibration products, e.g. combined bias, dark, flat, etc.", 

64 ) 

65 doAdd2DBias = pexConfig.Field( 

66 dtype=bool, 

67 default=True, 

68 doc="Add 2D bias residual frame to data.", 

69 ) 

70 doAddBrightDefects = pexConfig.Field( 

71 dtype=bool, 

72 default=True, 

73 doc="Add bright defects (bad column) to data.", 

74 ) 

75 brightDefectLevel = pexConfig.Field( 

76 dtype=float, 

77 default=30000.0, 

78 doc="Bright defect level (electron).", 

79 ) 

80 doAddBadParallelOverscanColumn = pexConfig.Field( 

81 dtype=bool, 

82 default=True, 

83 doc="Add a bad column to the parallel overscan.", 

84 ) 

85 badParallelOverscanColumnLevel = pexConfig.Field( 

86 dtype=float, 

87 default=300000., 

88 doc="Bright parallel overscan column level (electron). Should be above saturation.", 

89 ) 

90 doAddBadParallelOverscanColumnNeighbors = pexConfig.Field( 

91 dtype=bool, 

92 default=True, 

93 doc="Add low-level bad columns next to parallel overscan bad column.", 

94 ) 

95 badParallelOverscanColumnNeighborsLevel = pexConfig.Field( 

96 dtype=float, 

97 default=50.0, 

98 doc="Bright parallel overscan column neighbors level (electron).", 

99 ) 

100 doAddBrighterFatter = pexConfig.Field( 

101 dtype=bool, 

102 default=False, 

103 doc="Add brighter fatter and/or diffusion effects to image.", 

104 ) 

105 doAddDeferredCharge = pexConfig.Field( 

106 dtype=bool, 

107 default=False, 

108 doc="Add serial CTI at the amp level?", 

109 ) 

110 bfStrength = pexConfig.Field( 

111 dtype=float, 

112 default=2.0, 

113 doc="The brighter fatter effect scaling parameter (cannot be zero)." 

114 "Nominally = 1, but = 2 is more realistic." 

115 ) 

116 nRecalc = pexConfig.Field( 

117 dtype=int, 

118 default=10000, 

119 doc="Number of electrons to accumulate before recalculating pixel shapes.", 

120 ) 

121 doAddClockInjectedOffset = pexConfig.Field( 

122 dtype=bool, 

123 default=True, 

124 doc="Add clock-injected offset to data (on-chip bias level).", 

125 ) 

126 clockInjectedOffsetLevel = pexConfig.Field( 

127 dtype=float, 

128 default=8500.0, 

129 doc="Clock-injected offset (on-chip bias level), in electron.", 

130 ) 

131 noise2DBias = pexConfig.Field( 

132 dtype=float, 

133 default=2.0, 

134 doc="Noise (in electron) to generate a 2D bias residual frame.", 

135 ) 

136 doAddDarkNoiseOnly = pexConfig.Field( 

137 dtype=bool, 

138 default=False, 

139 doc="Add only dark current noise, for testing consistency.", 

140 ) 

141 doAddParallelOverscanRamp = pexConfig.Field( 

142 dtype=bool, 

143 default=True, 

144 doc="Add overscan ramp to parallel overscan and data regions.", 

145 ) 

146 doAddSerialOverscanRamp = pexConfig.Field( 

147 dtype=bool, 

148 default=True, 

149 doc="Add overscan ramp to serial overscan and data regions.", 

150 ) 

151 doAddHighSignalNonlinearity = pexConfig.Field( 

152 dtype=bool, 

153 default=True, 

154 doc="Add high signal non-linearity to overscan and data regions?", 

155 ) 

156 doAddLowSignalNonlinearity = pexConfig.Field( 

157 dtype=bool, 

158 default=False, 

159 doc="Add low signal non-linearity to overscan and data regions? (Not supported yet.)", 

160 ) 

161 highSignalNonlinearityThreshold = pexConfig.Field( 

162 dtype=float, 

163 default=40_000., 

164 doc="Threshold (in adu) for the non-linearity to be considered ``high signal``.", 

165 ) 

166 doApplyGain = pexConfig.Field( 

167 dtype=bool, 

168 default=True, 

169 doc="Add gain to data.", 

170 ) 

171 doRoundAdu = pexConfig.Field( 

172 dtype=bool, 

173 default=True, 

174 doc="Round adu values to nearest integer.", 

175 ) 

176 gainDict = pexConfig.DictField( 

177 keytype=str, 

178 itemtype=float, 

179 doc="Dictionary of amp name to gain; any amps not listed will use " 

180 "config.gain as the value. Units are electron/adu.", 

181 default={ 

182 "C:0,0": 1.65, 

183 "C:0,1": 1.60, 

184 "C:0,2": 1.55, 

185 "C:0,3": 1.70, 

186 "C:1,0": 1.75, 

187 "C:1,1": 1.80, 

188 "C:1,2": 1.85, 

189 "C:1,3": 1.70, 

190 }, 

191 ) 

192 assembleCcd = pexConfig.ConfigurableField( 

193 target=AssembleCcdTask, 

194 doc="CCD assembly task; used for defect box conversions.", 

195 ) 

196 

197 def validate(self): 

198 super().validate() 

199 

200 if self.doAddLowSignalNonlinearity: 

201 raise NotImplementedError("Low signal non-linearity is not implemented.") 

202 

203 if self.doAddDeferredCharge and self.isTrimmed: 

204 raise NotImplementedError("Must be untrimmed for mock serial CTI to" 

205 "realistically add charge into overscan regions.") 

206 

207 def setDefaults(self): 

208 super().setDefaults() 

209 

210 self.gain = 1.7 # Default value. 

211 self.skyLevel = 1700.0 # electron 

212 self.sourceFlux = [50_000.0] # electron 

213 self.sourceX = [35.0] # pixel 

214 self.sourceY = [37.0] # pixel 

215 self.overscanScale = 170.0 # electron 

216 self.biasLevel = 20_000.0 # adu 

217 self.doAddCrosstalk = True 

218 

219 

220class IsrMockLSST(IsrMock): 

221 """Class to generate consistent mock images for ISR testing. 

222 """ 

223 ConfigClass = IsrMockLSSTConfig 

224 _DefaultName = "isrMockLSST" 

225 

226 def __init__(self, **kwargs): 

227 super().__init__(**kwargs) 

228 

229 # Get kernel derived from imSim generated flats with BFE. The kernel 

230 # was used for Ops Rehearsal 3 for LSSTCam-type sensors 

231 # See https://rubinobs.atlassian.net/browse/DM-43059 for more details. 

232 self.bfKernel = np.array([[4.83499829e-01, 8.10171823e-01, 5.31096720e-01, 

233 3.54369868e-02, -8.44782871e-01, -1.64614462e+00, 

234 -3.83933101e+00, -5.60243416e+00, -6.51691578e+00, 

235 -5.60243416e+00, -3.83933101e+00, -1.64614462e+00, 

236 -8.44782871e-01, 3.54369868e-02, 5.31096720e-01, 

237 8.10171823e-01, 4.83499829e-01], 

238 [1.12382749e+00, 2.22609074e+00, 1.27877807e+00, 

239 4.55434098e-01, -1.76842385e+00, -1.90046460e+00, 

240 -8.10874526e+00, -1.20534899e+01, -1.48627948e+01, 

241 -1.20534899e+01, -8.10874526e+00, -1.90046460e+00, 

242 -1.76842385e+00, 4.55434098e-01, 1.27877807e+00, 

243 2.22609074e+00, 1.12382749e+00], 

244 [1.78571940e+00, 4.38918110e+00, 3.95098587e+00, 

245 3.70961649e-01, -3.48151981e+00, -9.61567736e+00, 

246 -1.78621172e+01, -2.32278872e+01, -2.31833727e+01, 

247 -2.32278872e+01, -1.78621172e+01, -9.61567736e+00, 

248 -3.48151981e+00, 3.70961649e-01, 3.95098587e+00, 

249 4.38918110e+00, 1.78571940e+00], 

250 [1.62986900e+00, 3.67851228e+00, 5.68645252e+00, 

251 2.15342566e-01, -8.89937202e+00, -1.44739813e+01, 

252 -2.98952660e+01, -4.37420817e+01, -4.83160958e+01, 

253 -4.37420817e+01, -2.98952660e+01, -1.44739813e+01, 

254 -8.89937202e+00, 2.15342566e-01, 5.68645252e+00, 

255 3.67851228e+00, 1.62986900e+00], 

256 [1.05524430e+00, 1.71917897e+00, 1.73105590e+00, 

257 -2.10088420e+00, -1.15118208e+01, -2.55007598e+01, 

258 -4.73056159e+01, -6.97257685e+01, -8.09264433e+01, 

259 -6.97257685e+01, -4.73056159e+01, -2.55007598e+01, 

260 -1.15118208e+01, -2.10088420e+00, 1.73105590e+00, 

261 1.71917897e+00, 1.05524430e+00], 

262 [8.71929228e-01, 5.41025574e-01, 9.47560771e-01, 

263 -5.75314708e-01, -7.46104027e+00, -4.42314481e+01, 

264 -9.54126971e+01, -1.61603201e+02, -2.07520692e+02, 

265 -1.61603201e+02, -9.54126971e+01, -4.42314481e+01, 

266 -7.46104027e+00, -5.75314708e-01, 9.47560771e-01, 

267 5.41025574e-01, 8.71929228e-01], 

268 [1.89144704e+00, 3.57543979e+00, -6.91419168e-02, 

269 -3.37950835e+00, -1.46695089e+01, -7.22850746e+01, 

270 -1.65563055e+02, -3.10820425e+02, -4.70026655e+02, 

271 -3.10820425e+02, -1.65563055e+02, -7.22850746e+01, 

272 -1.46695089e+01, -3.37950835e+00, -6.91419168e-02, 

273 3.57543979e+00, 1.89144704e+00], 

274 [3.11841913e+00, 7.84024994e+00, 1.88495248e+00, 

275 -7.69011009e+00, -2.71782400e+01, -1.04343326e+02, 

276 -2.47561370e+02, -5.32959841e+02, -1.16529012e+03, 

277 -5.32959841e+02, -2.47561370e+02, -1.04343326e+02, 

278 -2.71782400e+01, -7.69011009e+00, 1.88495248e+00, 

279 7.84024994e+00, 3.11841913e+00], 

280 [2.74197956e+00, 4.73107997e+00, -9.48352966e-01, 

281 -9.44822832e+00, -3.06477671e+01, -1.26788739e+02, 

282 -3.22828411e+02, -8.47943472e+02, -3.87702420e+03, 

283 -8.47943472e+02, -3.22828411e+02, -1.26788739e+02, 

284 -3.06477671e+01, -9.44822832e+00, -9.48352966e-01, 

285 4.73107997e+00, 2.74197956e+00], 

286 [3.11841913e+00, 7.84024994e+00, 1.88495248e+00, 

287 -7.69011009e+00, -2.71782400e+01, -1.04343326e+02, 

288 -2.47561370e+02, -5.32959841e+02, -1.16529012e+03, 

289 -5.32959841e+02, -2.47561370e+02, -1.04343326e+02, 

290 -2.71782400e+01, -7.69011009e+00, 1.88495248e+00, 

291 7.84024994e+00, 3.11841913e+00], 

292 [1.89144704e+00, 3.57543979e+00, -6.91419168e-02, 

293 -3.37950835e+00, -1.46695089e+01, -7.22850746e+01, 

294 -1.65563055e+02, -3.10820425e+02, -4.70026655e+02, 

295 -3.10820425e+02, -1.65563055e+02, -7.22850746e+01, 

296 -1.46695089e+01, -3.37950835e+00, -6.91419168e-02, 

297 3.57543979e+00, 1.89144704e+00], 

298 [8.71929228e-01, 5.41025574e-01, 9.47560771e-01, 

299 -5.75314708e-01, -7.46104027e+00, -4.42314481e+01, 

300 -9.54126971e+01, -1.61603201e+02, -2.07520692e+02, 

301 -1.61603201e+02, -9.54126971e+01, -4.42314481e+01, 

302 -7.46104027e+00, -5.75314708e-01, 9.47560771e-01, 

303 5.41025574e-01, 8.71929228e-01], 

304 [1.05524430e+00, 1.71917897e+00, 1.73105590e+00, 

305 -2.10088420e+00, -1.15118208e+01, -2.55007598e+01, 

306 -4.73056159e+01, -6.97257685e+01, -8.09264433e+01, 

307 -6.97257685e+01, -4.73056159e+01, -2.55007598e+01, 

308 -1.15118208e+01, -2.10088420e+00, 1.73105590e+00, 

309 1.71917897e+00, 1.05524430e+00], 

310 [1.62986900e+00, 3.67851228e+00, 5.68645252e+00, 

311 2.15342566e-01, -8.89937202e+00, -1.44739813e+01, 

312 -2.98952660e+01, -4.37420817e+01, -4.83160958e+01, 

313 -4.37420817e+01, -2.98952660e+01, -1.44739813e+01, 

314 -8.89937202e+00, 2.15342566e-01, 5.68645252e+00, 

315 3.67851228e+00, 1.62986900e+00], 

316 [1.78571940e+00, 4.38918110e+00, 3.95098587e+00, 

317 3.70961649e-01, -3.48151981e+00, -9.61567736e+00, 

318 -1.78621172e+01, -2.32278872e+01, -2.31833727e+01, 

319 -2.32278872e+01, -1.78621172e+01, -9.61567736e+00, 

320 -3.48151981e+00, 3.70961649e-01, 3.95098587e+00, 

321 4.38918110e+00, 1.78571940e+00], 

322 [1.12382749e+00, 2.22609074e+00, 1.27877807e+00, 

323 4.55434098e-01, -1.76842385e+00, -1.90046460e+00, 

324 -8.10874526e+00, -1.20534899e+01, -1.48627948e+01, 

325 -1.20534899e+01, -8.10874526e+00, -1.90046460e+00, 

326 -1.76842385e+00, 4.55434098e-01, 1.27877807e+00, 

327 2.22609074e+00, 1.12382749e+00], 

328 [4.83499829e-01, 8.10171823e-01, 5.31096720e-01, 

329 3.54369868e-02, -8.44782871e-01, -1.64614462+00, 

330 -3.83933101e+00, -5.60243416e+00, -6.51691578e+00, 

331 -5.60243416e+00, -3.83933101e+00, -1.64614462e+00, 

332 -8.44782871e-01, 3.54369868e-02, 5.31096720e-01, 

333 8.10171823e-01, 4.83499829e-01]]) * 1e-10 

334 

335 # Spline trap coefficients and the ctiCalibDict are all taken from a 

336 # cti calibration measured from LSSTCam sensor R03_S12 during Run 5 

337 # EO testing. These are the coefficients for the spline trap model 

338 # used in the deferred charge calibration. The collection can be 

339 # found in /repo/ir2: u/abrought/13144/cti (processed 3/4/2024). 

340 self.splineTrapCoeffs = np.array([0.0, 28.1, 47.4, 56.4, 66.6, 78.6, 92.4, 109.4, 

341 129.0, 151.9, 179.4, 211.9, 250.5, 296.2, 350.0, 

342 413.5, 488.0, 576.0, 680.4, 753.0, 888.2, 1040.5, 

343 1254.1, 1478.9, 1747.0, 2055.7, 2416.9, 2855.2, 

344 3361.9, 3969.4, 4665.9, 5405.3, 6380.0, 7516.7, 

345 8875.9, 10488.6, 12681.9, 14974.2, 17257.6, 20366.5, 

346 24026.7, 28372.1, 33451.7, 39550.4, 46624.8, 55042.9, 

347 64862.7, 76503.1, 90265.6, 106384.2, 0.0, 0.0, 0.0, 

348 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.1, 

349 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.1, 0.0, 

350 0.1, 0.2, 0.6, 0.1, 0.0, 0.1, 0.0, 0.6, 0.3, 0.5, 0.8, 

351 0.8, 1.5, 2.0, 1.8, 2.4, 2.6, 3.7, 5.0, 6.4, 8.4, 10.9, 

352 14.5, 21.1, 28.9]) 

353 

354 self.ctiCalibDict = {'driftScale': {'C:0,0': 0.000127, 

355 'C:0,1': 0.000137, 

356 'C:0,2': 0.000138, 

357 'C:0,3': 0.000147, 

358 'C:1,0': 0.000147, 

359 'C:1,1': 0.000122, 

360 'C:1,2': 0.000123, 

361 'C:1,3': 0.000116}, 

362 'decayTime': {'C:0,0': 2.30, 

363 'C:0,1': 2.21, 

364 'C:0,2': 2.28, 

365 'C:0,3': 2.34, 

366 'C:1,0': 2.30, 

367 'C:1,1': 2.40, 

368 'C:1,2': 2.51, 

369 'C:1,3': 2.21}, 

370 'globalCti': {'C:0,0': 5.25e-07, 

371 'C:0,1': 5.38e-07, 

372 'C:0,2': 5.80e-07, 

373 'C:0,3': 5.91e-07, 

374 'C:1,0': 6.24e-07, 

375 'C:1,1': 5.72e-07, 

376 'C:1,2': 5.60e-07, 

377 'C:1,3': 4.40e-07}, 

378 'signals': {'C:0,0': np.linspace(1.0e2, 1.0e5, 10), 

379 'C:0,1': np.linspace(1.0e2, 1.0e5, 10), 

380 'C:0,2': np.linspace(1.0e2, 1.0e5, 10), 

381 'C:0,3': np.linspace(1.0e2, 1.0e5, 10), 

382 'C:1,0': np.linspace(1.0e2, 1.0e5, 10), 

383 'C:1,1': np.linspace(1.0e2, 1.0e5, 10), 

384 'C:1,2': np.linspace(1.0e2, 1.0e5, 10), 

385 'C:1,3': np.linspace(1.0e2, 1.0e5, 10)}, 

386 'serialEper': {'C:0,0': np.full(10, 5.25e-07), 

387 'C:0,1': np.full(10, 5.38e-07), 

388 'C:0,2': np.full(10, 5.80e-07), 

389 'C:0,3': np.full(10, 5.91e-07), 

390 'C:1,0': np.full(10, 6.24e-07), 

391 'C:1,1': np.full(10, 5.72e-07), 

392 'C:1,2': np.full(10, 5.60e-07), 

393 'C:1,3': np.full(10, 4.40e-07)}, 

394 'parallelEper': {'C:0,0': np.full(10, 5.25e-07), 

395 'C:0,1': np.full(10, 5.38e-07), 

396 'C:0,2': np.full(10, 5.80e-07), 

397 'C:0,3': np.full(10, 5.91e-07), 

398 'C:1,0': np.full(10, 6.24e-07), 

399 'C:1,1': np.full(10, 5.72e-07), 

400 'C:1,2': np.full(10, 5.60e-07), 

401 'C:1,3': np.full(10, 4.40e-07)}, 

402 'serialCtiTurnoff': {'C:0,0': 1.0e5, 

403 'C:0,1': 1.0e5, 

404 'C:0,2': 1.0e5, 

405 'C:0,3': 1.0e5, 

406 'C:1,0': 1.0e5, 

407 'C:1,1': 1.0e5, 

408 'C:1,2': 1.0e5, 

409 'C:1,3': 1.0e5}, 

410 'parallelCtiTurnoff': {'C:0,0': 1.0e5, 

411 'C:0,1': 1.0e5, 

412 'C:0,2': 1.0e5, 

413 'C:0,3': 1.0e5, 

414 'C:1,0': 1.0e5, 

415 'C:1,1': 1.0e5, 

416 'C:1,2': 1.0e5, 

417 'C:1,3': 1.0e5}, 

418 'serialCtiTurnoffSamplingErr': {'C:0,0': 1.0e3, 

419 'C:0,1': 1.0e3, 

420 'C:0,2': 1.0e3, 

421 'C:0,3': 1.0e3, 

422 'C:1,0': 1.0e3, 

423 'C:1,1': 1.0e3, 

424 'C:1,2': 1.0e3, 

425 'C:1,3': 1.0e3}, 

426 'parallelCtiTurnoffSamplingErr': {'C:0,0': 1.0e3, 

427 'C:0,1': 1.0e3, 

428 'C:0,2': 1.0e3, 

429 'C:0,3': 1.0e3, 

430 'C:1,0': 1.0e3, 

431 'C:1,1': 1.0e3, 

432 'C:1,2': 1.0e3, 

433 'C:1,3': 1.0e3}, 

434 'serialTraps': {'C:0,0': {'size': 20000.0, 

435 'emissionTime': 0.4, 

436 'pixel': 1, 

437 'trap_type': 'spline', 

438 'coeffs': self.splineTrapCoeffs}, 

439 'C:0,1': {'size': 20000.0, 

440 'emissionTime': 0.4, 

441 'pixel': 1, 

442 'trap_type': 'spline', 

443 'coeffs': self.splineTrapCoeffs}, 

444 'C:0,2': {'size': 20000.0, 

445 'emissionTime': 0.4, 

446 'pixel': 1, 

447 'trap_type': 'spline', 

448 'coeffs': self.splineTrapCoeffs}, 

449 'C:0,3': {'size': 20000.0, 

450 'emissionTime': 0.4, 

451 'pixel': 1, 

452 'trap_type': 'spline', 

453 'coeffs': self.splineTrapCoeffs}, 

454 'C:1,0': {'size': 20000.0, 

455 'emissionTime': 0.4, 

456 'pixel': 1, 

457 'trap_type': 'spline', 

458 'coeffs': self.splineTrapCoeffs}, 

459 'C:1,1': {'size': 20000.0, 

460 'emissionTime': 0.4, 

461 'pixel': 1, 

462 'trap_type': 'spline', 

463 'coeffs': self.splineTrapCoeffs}, 

464 'C:1,2': {'size': 20000.0, 

465 'emissionTime': 0.4, 

466 'pixel': 1, 

467 'trap_type': 'spline', 

468 'coeffs': self.splineTrapCoeffs}, 

469 'C:1,3': {'size': 20000.0, 

470 'emissionTime': 0.4, 

471 'pixel': 1, 

472 'trap_type': 'spline', 

473 'coeffs': self.splineTrapCoeffs}}} 

474 

475 self.deferredChargeCalib = self.makeDeferredChargeCalib() 

476 

477 # Cross-talk coeffs are defined in the parent class. 

478 

479 self.makeSubtask("assembleCcd") 

480 

481 def run(self): 

482 """Generate a mock ISR product following LSSTCam ISR, and return it. 

483 

484 Returns 

485 ------- 

486 image : `lsst.afw.image.Exposure` 

487 Simulated ISR image with signals added. 

488 dataProduct : 

489 Simulated ISR data products. 

490 None : 

491 Returned if no valid configuration was found. 

492 

493 Raises 

494 ------ 

495 RuntimeError 

496 Raised if both doGenerateImage and doGenerateData are specified. 

497 """ 

498 if self.config.doGenerateImage and self.config.doGenerateData: 

499 raise RuntimeError("Only one of doGenerateImage and doGenerateData may be specified.") 

500 elif self.config.doGenerateImage: 

501 return self.makeImage() 

502 elif self.config.doGenerateData: 

503 return self.makeData() 

504 else: 

505 return None 

506 

507 def makeImage(self): 

508 """Generate a simulated ISR LSST image. 

509 

510 Returns 

511 ------- 

512 exposure : `lsst.afw.image.Exposure` or `dict` 

513 Simulated ISR image data. 

514 

515 Notes 

516 ----- 

517 This method constructs a "raw" data image. 

518 """ 

519 exposure = self.getExposure() 

520 

521 # Set up random number generators for consistency of components, 

522 # no matter the group that are configured. 

523 rngSky = np.random.RandomState(seed=self.config.rngSeed + 1) 

524 rngDark = np.random.RandomState(seed=self.config.rngSeed + 2) 

525 rng2DBias = np.random.RandomState(seed=self.config.rngSeed + 3) 

526 rngOverscan = np.random.RandomState(seed=self.config.rngSeed + 4) 

527 rngReadNoise = np.random.RandomState(seed=self.config.rngSeed + 5) 

528 rngBrighterFatter = galsim.BaseDeviate(self.config.rngSeed + 6) 

529 

530 # Create the linearizer if we will need it. 

531 if self.config.doAddHighSignalNonlinearity: 

532 linearizer = LinearizerMockLSST().run() 

533 

534 # We introduce effects as they happen from a source to the signal, 

535 # so the effects go from electron to adu. 

536 # The ISR steps will then correct these effects in the reverse order. 

537 for idx, amp in enumerate(exposure.getDetector()): 

538 # Get image bbox and data 

539 bbox = None 

540 if self.config.isTrimmed: 

541 bbox = amp.getBBox() 

542 bboxFull = bbox 

543 else: 

544 bbox = amp.getRawDataBBox() 

545 bboxFull = amp.getRawBBox() 

546 

547 # This is the image data (excluding pre/overscans). 

548 ampImageData = exposure.image[bbox] 

549 # This is the full data (including pre/overscans if untrimmed). 

550 ampFullData = exposure.image[bboxFull] 

551 

552 # Astrophysical signals are all in electron (e-). 

553 # These are only applied to the imaging portion of the 

554 # amplifier (ampImageData) 

555 

556 if self.config.doAddSky: 

557 # The sky effects are in electron. 

558 self.amplifierAddNoise( 

559 ampImageData, 

560 self.config.skyLevel, 

561 np.sqrt(self.config.skyLevel), 

562 rng=rngSky, 

563 ) 

564 

565 if self.config.doAddSource: 

566 for sourceAmp, sourceFlux, sourceX, sourceY in zip(self.config.sourceAmp, 

567 self.config.sourceFlux, 

568 self.config.sourceX, 

569 self.config.sourceY): 

570 if idx == sourceAmp: 

571 # The source flux is in electron. 

572 self.amplifierAddSource(ampImageData, sourceFlux, sourceX, sourceY) 

573 

574 if self.config.doAddFringe: 

575 # Fringes are added in electron. 

576 self.amplifierAddFringe(amp, 

577 ampImageData, 

578 np.array(self.config.fringeScale), 

579 x0=np.array(self.config.fringeX0), 

580 y0=np.array(self.config.fringeY0)) 

581 

582 if self.config.doAddFlat: 

583 if self.config.calibMode: 

584 # In case we are making a combined flat, 

585 # add a non-zero signal so the mock flat can be multiplied 

586 self.amplifierAddNoise(ampImageData, 1.0, 0.0) 

587 # Multiply each amplifier by a Gaussian centered on u0 and v0 

588 u0 = exposure.getDetector().getBBox().getDimensions().getX()/2. 

589 v0 = exposure.getDetector().getBBox().getDimensions().getY()/2. 

590 self.amplifierMultiplyFlat(amp, ampImageData, self.config.flatDrop, u0=u0, v0=v0) 

591 

592 # On-chip electronic effects. 

593 

594 # 1. Add bright defect(s). 

595 if self.config.doAddBrightDefects: 

596 defectList = self.makeDefectList(isTrimmed=self.config.isTrimmed) 

597 

598 for defect in defectList: 

599 exposure.image[defect.getBBox()] = self.config.brightDefectLevel 

600 

601 for idx, amp in enumerate(exposure.getDetector()): 

602 # Get image bbox and data 

603 bbox = None 

604 if self.config.isTrimmed: 

605 bbox = amp.getBBox() 

606 bboxFull = bbox 

607 else: 

608 bbox = amp.getRawDataBBox() 

609 bboxFull = amp.getRawBBox() 

610 

611 # This is the image data (excluding pre/overscans). 

612 ampImageData = exposure.image[bbox] 

613 # This is the full data (including pre/overscans if untrimmed). 

614 ampFullData = exposure.image[bboxFull] 

615 

616 # 2. Add dark current (electron) to imaging portion of the amp. 

617 if self.config.doAddDark or self.config.doAddDarkNoiseOnly: 

618 if self.config.doAddDarkNoiseOnly: 

619 darkLevel = 0.0 

620 else: 

621 darkLevel = self.config.darkRate * self.config.darkTime 

622 if self.config.calibMode: 

623 darkNoise = 0.0 

624 else: 

625 darkNoise = np.sqrt(self.config.darkRate * self.config.darkTime) 

626 

627 self.amplifierAddNoise(ampImageData, darkLevel, darkNoise, rng=rngDark) 

628 

629 # 3. Add BF effect (electron) to imaging portion of the amp. 

630 if self.config.doAddBrighterFatter is True: 

631 self.amplifierAddBrighterFatter( 

632 ampImageData, 

633 rngBrighterFatter, 

634 self.config.bfStrength, 

635 self.config.nRecalc, 

636 ) 

637 

638 # 4. Add serial CTI (electron) to amplifier (imaging + overscan). 

639 if self.config.doAddDeferredCharge: 

640 # Get the free charge area for the amplifier. 

641 self.amplifierAddDeferredCharge(exposure, amp) 

642 

643 # 5. Add 2D bias residual (electron) to imaging portion of the amp. 

644 if self.config.doAdd2DBias: 

645 # For now we use an unstructured noise field to add some 

646 # consistent 2D bias residual that can be subtracted. In 

647 # the future this can be made into a warm corner (for example). 

648 self.amplifierAddNoise( 

649 ampImageData, 

650 0.0, 

651 self.config.noise2DBias, 

652 rng=rng2DBias, 

653 ) 

654 

655 # 6. Add clock-injected offset (electron) to amplifer 

656 # (imaging + overscan). 

657 # This is just an offset that will be crosstalked and modified by 

658 # the gain, and does not have a noise associated with it. 

659 if self.config.doAddClockInjectedOffset: 

660 self.amplifierAddNoise( 

661 ampFullData, 

662 self.config.clockInjectedOffsetLevel, 

663 0.0, 

664 ) 

665 

666 # 7./8. Add serial and parallel overscan slopes (electron) 

667 # (imaging + overscan) 

668 if (self.config.doAddParallelOverscanRamp or self.config.doAddSerialOverscanRamp) and \ 

669 not self.config.isTrimmed: 

670 

671 if self.config.doAddParallelOverscanRamp: 

672 # Apply gradient along the X axis. 

673 self.amplifierAddXGradient(ampFullData, -1.0 * self.config.overscanScale, 

674 1.0 * self.config.overscanScale) 

675 

676 if self.config.doAddSerialOverscanRamp: 

677 # Apply the gradient along the Y axis. 

678 self.amplifierAddYGradient(ampFullData, -1.0 * self.config.overscanScale, 

679 1.0 * self.config.overscanScale) 

680 

681 # 9. Add non-linearity (electron) to amplifier 

682 # (imaging + overscan). 

683 if self.config.doAddHighSignalNonlinearity: 

684 # The linearizer coefficients come from makeLinearizer(). 

685 if linearizer.linearityType[amp.getName()] != "Spline": 

686 raise RuntimeError("IsrMockLSST only supports spline non-linearity.") 

687 

688 coeffs = linearizer.linearityCoeffs[amp.getName()] 

689 centers, values = np.split(coeffs, 2) 

690 

691 # This is an application of high signal non-linearity, so we 

692 # set the lower values to 0.0 (this cut is arbitrary). 

693 values[centers < self.config.highSignalNonlinearityThreshold] = 0.0 

694 

695 # The linearizer is units of adu, so convert to electron 

696 centers *= self.config.gainDict[amp.getName()] 

697 values *= self.config.gainDict[amp.getName()] 

698 

699 # Note that the linearity spline is in "overscan subtracted" 

700 # units so needs to be applied without the clock-injected 

701 # offset. 

702 self.amplifierAddNonlinearity( 

703 ampFullData, 

704 centers, 

705 values, 

706 self.config.clockInjectedOffsetLevel if self.config.doAddClockInjectedOffset else 0.0, 

707 ) 

708 

709 # 10. Add read noise (electron) to the amplifier 

710 # (imaging + overscan). 

711 # Unsure if this should be before or after crosstalk. 

712 # Probably some of both; hopefully doesn't matter. 

713 if not self.config.calibMode: 

714 # Add read noise to the imaging region. 

715 self.amplifierAddNoise( 

716 ampImageData, 

717 0.0, 

718 self.config.readNoise, 

719 rng=rngReadNoise, 

720 ) 

721 

722 # If not trimmed, add to the overscan regions. 

723 if not self.config.isTrimmed: 

724 parallelOverscanBBox = amp.getRawParallelOverscanBBox() 

725 parallelOverscanData = exposure.image[parallelOverscanBBox] 

726 

727 serialOverscanBBox = self.getFullSerialOverscanBBox(amp) 

728 serialOverscanData = exposure.image[serialOverscanBBox] 

729 

730 # Add read noise of mean 0 

731 # to the parallel and serial overscan regions. 

732 self.amplifierAddNoise( 

733 parallelOverscanData, 

734 0.0, 

735 self.config.readNoise, 

736 rng=rngOverscan, 

737 ) 

738 self.amplifierAddNoise( 

739 serialOverscanData, 

740 0.0, 

741 self.config.readNoise, 

742 rng=rngOverscan, 

743 ) 

744 

745 # 7b. Add bad column to the parallel overscan region. 

746 if self.config.doAddBadParallelOverscanColumn and not self.config.isTrimmed: 

747 # We want to place this right above the defect, to simulate 

748 # bleeding into the parallel overscan region. 

749 amp = exposure.getDetector()[2] 

750 parBBox = amp.getRawParallelOverscanBBox() 

751 bboxBad = geom.Box2I( 

752 corner=geom.Point2I(50, parBBox.getMinY()), 

753 dimensions=geom.Extent2I(1, parBBox.getHeight()), 

754 ) 

755 exposure[bboxBad].image.array[:, :] += self.config.badParallelOverscanColumnLevel 

756 

757 if self.config.doAddBadParallelOverscanColumnNeighbors: 

758 for neighbor in [49, 51]: 

759 bboxBad = geom.Box2I( 

760 corner=geom.Point2I(neighbor, parBBox.getMinY()), 

761 dimensions=geom.Extent2I(1, parBBox.getHeight()), 

762 ) 

763 exposure[bboxBad].image.array[:, :] += self.config.badParallelOverscanColumnNeighborsLevel 

764 

765 # 11. Add crosstalk (electron) to all the amplifiers 

766 # (imaging + overscan). 

767 if self.config.doAddCrosstalk: 

768 ctCalib = CrosstalkCalib() 

769 # We use the regular subtractCrosstalk code but with a negative 

770 # sign on the crosstalk coefficients so it adds instead of 

771 # subtracts. We only apply the signal plane (ignoreVariance, 

772 # subtrahendMasking) with a very large pixel to mask to ensure 

773 # no crosstalk mask bits are set. 

774 ctCalib.subtractCrosstalk( 

775 exposure, 

776 crosstalkCoeffs=-1*self.crosstalkCoeffs, 

777 doSubtrahendMasking=True, 

778 minPixelToMask=np.inf, 

779 ignoreVariance=True, 

780 fullAmplifier=True, 

781 ) 

782 

783 for amp in exposure.getDetector(): 

784 # Get image bbox and data (again). 

785 bbox = None 

786 if self.config.isTrimmed: 

787 bbox = amp.getBBox() 

788 bboxFull = bbox 

789 else: 

790 bbox = amp.getRawDataBBox() 

791 bboxFull = amp.getRawBBox() 

792 

793 # This is the image data (excluding pre/overscans). 

794 ampImageData = exposure.image[bbox] 

795 # This is the full data (including pre/overscans if untrimmed). 

796 ampFullData = exposure.image[bboxFull] 

797 

798 # 12. Gain un-normalize (from electron to floating point adu) 

799 if self.config.doApplyGain: 

800 gain = self.config.gainDict.get(amp.getName(), self.config.gain) 

801 self.applyGain(ampFullData, gain) 

802 

803 # 13. Add overall bias level (adu) to the amplifier 

804 # (imaging + overscan) 

805 if self.config.doAddBias: 

806 self.addBiasLevel(ampFullData, self.config.biasLevel) 

807 

808 # 14. Round/Truncate to integers (adu) 

809 if self.config.doRoundAdu: 

810 self.roundADU(ampFullData) 

811 

812 # Add units metadata to calibrations. 

813 if self.config.calibMode: 

814 if self.config.doApplyGain: 

815 exposure.metadata["LSST ISR UNITS"] = "adu" 

816 else: 

817 exposure.metadata["LSST ISR UNITS"] = "electron" 

818 

819 # Add a variance plane appropriate for a calibration frame. 

820 # We take the absolute value for biases which have no signal. 

821 exposure.variance.array[:, :] = np.abs(np.median(exposure.image.array)/10.) 

822 

823 exposure.metadata["BSSVBS"] = 50.0 

824 exposure.metadata["HVBIAS"] = "ON" 

825 

826 if self.config.doGenerateAmpDict: 

827 expDict = dict() 

828 for amp in exposure.getDetector(): 

829 expDict[amp.getName()] = exposure 

830 return expDict 

831 else: 

832 return exposure 

833 

834 def addBiasLevel(self, ampData, biasLevel): 

835 """Add bias level to an amplifier's image data. 

836 

837 Parameters 

838 ---------- 

839 ampData : `lsst.afw.image.ImageF` 

840 Amplifier image to operate on. 

841 biasLevel : `float` 

842 Bias level to be added to the image. 

843 """ 

844 ampArr = ampData.array 

845 ampArr[:] = ampArr[:] + biasLevel 

846 

847 def makeDefectList(self, isTrimmed=True): 

848 """Generate a simple defect list. 

849 

850 Parameters 

851 ---------- 

852 isTrimmed : `bool`, optional 

853 Return defects in trimmed coordinates? 

854 

855 Returns 

856 ------- 

857 defectList : `lsst.meas.algorithms.Defects` 

858 Simulated defect list 

859 """ 

860 defectBoxesUntrimmed = [ 

861 geom.Box2I( 

862 geom.Point2I(50, 118), 

863 geom.Extent2I(1, 51), 

864 ), 

865 ] 

866 

867 if not isTrimmed: 

868 return Defects(defectBoxesUntrimmed) 

869 

870 # If trimmed, we need to convert. 

871 tempExp = self.getExposure(isTrimmed=False) 

872 tempExp.image.array[:, :] = 0.0 

873 for bbox in defectBoxesUntrimmed: 

874 tempExp.image[bbox] = 1.0 

875 

876 assembledExp = self.assembleCcd.assembleCcd(tempExp) 

877 

878 # Use thresholding code to find defect footprints/boxes. 

879 threshold = afwDetection.createThreshold(1.0, "value", polarity=True) 

880 footprintSet = afwDetection.FootprintSet(assembledExp.image, threshold) 

881 

882 return Defects.fromFootprintList(footprintSet.getFootprints()) 

883 

884 def makeBfKernel(self): 

885 """Generate a simple simulated brighter-fatter kernel. 

886 Returns 

887 ------- 

888 kernel : `lsst.ip.isr.BrighterFatterKernel` 

889 Simulated brighter-fatter kernel. 

890 """ 

891 bfkArray = super().makeBfKernel() 

892 bfKernelObject = BrighterFatterKernel() 

893 bfKernelObject.level = 'AMP' 

894 bfKernelObject.gain = self.config.gainDict 

895 

896 for amp in self.getExposure().getDetector(): 

897 # Kernel must be in (y,x) orientation 

898 bfKernelObject.ampKernels[amp.getName()] = bfkArray.T 

899 bfKernelObject.valid[amp.getName()] = True 

900 

901 return bfKernelObject 

902 

903 def makeDeferredChargeCalib(self): 

904 """Generate a CTI calibration. 

905 

906 Returns 

907 ------- 

908 cti : `lsst.ip.isr.deferredCharge.DeferredChargeCalib` 

909 Simulated deferred charge calibration. 

910 """ 

911 

912 metadataDict = {'metadata': PropertyList()} 

913 metadataDict['metadata'].add(name="OBSTYPE", value="CTI") 

914 metadataDict['metadata'].add(name="CALIBCLS", 

915 value="lsst.ip.isr.deferredCharge.DeferredChargeCalib") 

916 self.ctiCalibDict = {**metadataDict, **self.ctiCalibDict} 

917 deferredChargeCalib = DeferredChargeCalib() 

918 self.cti = deferredChargeCalib.fromDict(self.ctiCalibDict) 

919 

920 return self.cti 

921 

922 def amplifierAddBrighterFatter(self, ampImageData, rng, bfStrength, nRecalc): 

923 """Add brighter fatter effect and/or diffusion to the image. 

924 Parameters 

925 ---------- 

926 ampImageData : `lsst.afw.image.ImageF` 

927 Amplifier image to operate on. 

928 rng : `galsim.BaseDeviate` 

929 Random number generator. 

930 bfStrength : `float` 

931 Scaling parameter of the brighter fatter effect (nominally = 1) 

932 nRecalc: 'int' 

933 The number of electrons to accumulate before recalculating the 

934 distortion of the pixel shapes. 

935 """ 

936 

937 incidentImage = galsim.Image(ampImageData.array, scale=1) 

938 measuredImage = galsim.ImageF( 

939 ampImageData.array.shape[1], 

940 ampImageData.array.shape[0], 

941 scale=1, 

942 ) 

943 photons = galsim.PhotonArray.makeFromImage(incidentImage) 

944 

945 sensorModel = galsim.SiliconSensor( 

946 strength=bfStrength, 

947 rng=rng, 

948 diffusion_factor=0.0, 

949 nrecalc=nRecalc, 

950 ) 

951 

952 totalFluxAdded = sensorModel.accumulate(photons, measuredImage) 

953 ampImageData.array = measuredImage.array 

954 

955 return totalFluxAdded 

956 

957 def amplifierAddDeferredCharge(self, exposure, amp): 

958 """Add serial CTI to the amplifier data. 

959 

960 Parameters 

961 ---------- 

962 exposure : `lsst.afw.image.ExposureF` 

963 The exposure object containing the amplifier 

964 to apply deferred charge to. 

965 amp : `lsst.afw.image.Amplifier` 

966 The amplifier object (contains geometry info). 

967 """ 

968 # Get the amplifier's geometry parameters. 

969 # When adding deferred charge, we have already assured that 

970 # isTrimmed is False. Therefore we want to make sure that we 

971 # get the RawDataBBox. 

972 readoutCorner = amp.getReadoutCorner() 

973 prescanWidth = amp.getRawHorizontalPrescanBBox().getWidth() 

974 serialOverscanWidth = amp.getRawHorizontalOverscanBBox().getWidth() 

975 parallelOverscanWidth = amp.getRawVerticalOverscanBBox().getHeight() 

976 bboxFreeCharge = amp.getRawDataBBox() 

977 bboxFreeCharge = bboxFreeCharge.expandedTo(amp.getRawHorizontalOverscanBBox()) 

978 bboxFreeCharge = bboxFreeCharge.expandedTo(amp.getRawHorizontalPrescanBBox()) 

979 bboxFreeCharge = bboxFreeCharge.expandedTo(amp.getRawVerticalOverscanBBox()) 

980 

981 ampFreeChargeData = exposure.image[bboxFreeCharge] 

982 ampImageData = exposure.image[amp.getRawDataBBox()] 

983 

984 # Get the deferred charge parameters for this amplifier. 

985 cti = self.deferredChargeCalib.globalCti[amp.getName()] 

986 traps = self.deferredChargeCalib.serialTraps[amp.getName()] 

987 driftScale = self.deferredChargeCalib.driftScale[amp.getName()] 

988 decayTime = self.deferredChargeCalib.decayTime[amp.getName()] 

989 

990 # Create a fake amplifier object that contains some deferred charge 

991 # paramters. 

992 floatingOutputAmplifier = FloatingOutputAmplifier( 

993 gain=1.0, # Everything is already in electrons. 

994 scale=driftScale, 

995 decay_time=decayTime, 

996 noise=0.0, 

997 offset=0.0, 

998 ) 

999 

1000 def flipImage(arr, readoutCorner): 

1001 # Flip an array so that the readout corner is in 

1002 # the lower left. 

1003 if readoutCorner == ReadoutCorner.LR: 

1004 return np.fliplr(arr) 

1005 elif readoutCorner == ReadoutCorner.UR: 

1006 return np.fliplr(np.flipud(arr)) 

1007 elif readoutCorner == ReadoutCorner.UL: 

1008 return np.flipud(arr) 

1009 else: 

1010 pass 

1011 

1012 return arr 

1013 

1014 # The algorithm expects that the readout corner is in 

1015 # the lower left corner. Flip it to be so: 

1016 imageData = ampImageData.array 

1017 imageData = flipImage(imageData, readoutCorner) 

1018 

1019 # Simulate the amplifier. 

1020 ampSim = SegmentSimulator( 

1021 imarr=imageData, 

1022 prescan_width=prescanWidth, 

1023 output_amplifier=floatingOutputAmplifier, 

1024 cti=cti, 

1025 traps=traps, 

1026 ) 

1027 

1028 # Simulate deferred charge! 

1029 # Note that the readout() method uses the image region data and the 

1030 # overscan dimensions provided as input parameters. It then creates 

1031 # overscan nd adds it to the image data to create the raw image. 

1032 result = ampSim.readout( 

1033 serial_overscan_width=serialOverscanWidth, 

1034 parallel_overscan_width=parallelOverscanWidth, 

1035 ) 

1036 

1037 # Flip the image back to the original orientation. 

1038 result = flipImage(result, readoutCorner) 

1039 

1040 # Set the image with the deferred charge added. 

1041 ampFreeChargeData.array[:, :] = result 

1042 

1043 def makeLinearizer(self): 

1044 # docstring inherited. 

1045 

1046 # The linearizer has units of adu. 

1047 nNodes = 10 

1048 # Set this to just above the mock saturation (adu) 

1049 maxADU = 101_000 

1050 nonLinSplineNodes = np.linspace(0, maxADU, nNodes) 

1051 # These values come from cp_pipe/tests/test_linearity.py and 

1052 # are based on a test fit to LSSTCam data, run 7193D, detector 22, 

1053 # amp C00. 

1054 nonLinSplineValues = np.array( 

1055 [0.0, -8.87, 1.46, 1.69, -6.92, -68.23, -78.01, -11.56, 80.26, 185.01] 

1056 ) 

1057 

1058 if self.config.doAddHighSignalNonlinearity and not self.config.doAddLowSignalNonlinearity: 

1059 nonLinSplineValues[nonLinSplineNodes < self.config.highSignalNonlinearityThreshold] = 0.0 

1060 elif self.config.doAddLowSignalNonlinearity: 

1061 raise NotImplementedError("Low signal non-linearity is not implemented.") 

1062 

1063 exp = self.getExposure() 

1064 detector = exp.getDetector() 

1065 

1066 linearizer = Linearizer(detector=detector) 

1067 linearizer.updateMetadataFromExposures([exp]) 

1068 

1069 # We need to set override by hand because we are constructing a 

1070 # linearizer manually and not from a serialized object. 

1071 linearizer.override = True 

1072 linearizer.hasLinearity = True 

1073 linearizer.validate() 

1074 linearizer.updateMetadata(camera=self.getCamera(), detector=detector, filterName='NONE') 

1075 linearizer.updateMetadata(setDate=True, setCalibId=True) 

1076 

1077 for amp in detector: 

1078 ampName = amp.getName() 

1079 linearizer.linearityType[ampName] = "Spline" 

1080 linearizer.linearityCoeffs[ampName] = np.concatenate([nonLinSplineNodes, nonLinSplineValues]) 

1081 # We need to specify the raw bbox here. 

1082 linearizer.linearityBBox[ampName] = amp.getRawBBox() 

1083 

1084 return linearizer 

1085 

1086 def amplifierAddNonlinearity(self, ampData, centers, values, offset): 

1087 """Add non-linearity to amplifier data. 

1088 

1089 Parameters 

1090 ---------- 

1091 ampData : `lsst.afw.image.ImageF` 

1092 Amplifier image to operate on. 

1093 centers : `np.ndarray` 

1094 Spline nodes. 

1095 values : `np.ndarray` 

1096 Spline values. 

1097 offset : `float` 

1098 Offset zero-point between linearizer (internal vs external). 

1099 """ 

1100 # I'm not sure what to do about negative values... 

1101 

1102 # Note that we are using the afw AKIMA_SPLINE to offset the 

1103 # data but using the equivalent but faster scipy Akima1DInterpolator to 

1104 # correct the data. 

1105 spl = afwMath.makeInterpolate( 

1106 centers, 

1107 values, 

1108 afwMath.stringToInterpStyle("AKIMA_SPLINE"), 

1109 ) 

1110 

1111 delta = np.asarray(spl.interpolate(ampData.array.ravel() - offset)) 

1112 

1113 ampData.array[:, :] += delta.reshape(ampData.array.shape) 

1114 

1115 def amplifierMultiplyFlat(self, amp, ampData, fracDrop, u0=100.0, v0=100.0): 

1116 """Multiply an amplifier's image data by a flat-like pattern. 

1117 

1118 Parameters 

1119 ---------- 

1120 amp : `lsst.afw.ampInfo.AmpInfoRecord` 

1121 Amplifier to operate on. Needed for amp<->exp coordinate 

1122 transforms. 

1123 ampData : `lsst.afw.image.ImageF` 

1124 Amplifier image to operate on. 

1125 fracDrop : `float` 

1126 Fractional drop from center to edge of detector along x-axis. 

1127 u0 : `float` 

1128 Peak location in detector coordinates. 

1129 v0 : `float` 

1130 Peak location in detector coordinates. 

1131 """ 

1132 if fracDrop >= 1.0: 

1133 raise RuntimeError("Flat fractional drop cannot be greater than 1.0") 

1134 

1135 sigma = u0 / np.sqrt(2.0 * fracDrop) 

1136 

1137 for x in range(0, ampData.getDimensions().getX()): 

1138 for y in range(0, ampData.getDimensions().getY()): 

1139 (u, v) = self.localCoordToExpCoord(amp, x, y) 

1140 f = np.exp(-0.5 * ((u - u0)**2 + (v - v0)**2) / sigma**2) 

1141 ampData.array[y][x] = (ampData.array[y][x] * f) 

1142 

1143 def applyGain(self, ampData, gain): 

1144 """Apply gain to the amplifier's data. 

1145 This method divides the data by the gain 

1146 because the mocks need to convert the data in electron to adu, 

1147 so it does the inverse operation to applyGains in isrFunctions. 

1148 

1149 Parameters 

1150 ---------- 

1151 ampData : `lsst.afw.image.ImageF` 

1152 Amplifier image to operate on. 

1153 gain : `float` 

1154 Gain value in electron/adu. 

1155 """ 

1156 ampArr = ampData.array 

1157 ampArr[:] = ampArr[:] / gain 

1158 

1159 def roundADU(self, ampData): 

1160 """Round adu to nearest integer. 

1161 

1162 Parameters 

1163 ---------- 

1164 ampData : `lsst.afw.image.ImageF` 

1165 Amplifier image to operate on. 

1166 """ 

1167 ampArr = ampData.array 

1168 ampArr[:] = np.around(ampArr) 

1169 

1170 def amplifierAddXGradient(self, ampData, start, end): 

1171 """Add a x-axis linear gradient to an amplifier's image data. 

1172 

1173 This method operates in the amplifier coordinate frame. 

1174 

1175 Parameters 

1176 ---------- 

1177 ampData : `lsst.afw.image.ImageF` 

1178 Amplifier image to operate on. 

1179 start : `float` 

1180 Start value of the gradient (at x=0). 

1181 end : `float` 

1182 End value of the gradient (at x=xmax). 

1183 """ 

1184 nPixX = ampData.getDimensions().getX() 

1185 ampArr = ampData.array 

1186 ampArr[:] = ampArr[:] + (np.interp(range(nPixX), (0, nPixX - 1), (start, end)).reshape(1, nPixX) 

1187 + np.zeros(ampData.getDimensions()).transpose()) 

1188 

1189 def getFullSerialOverscanBBox(self, amp): 

1190 """Get the full serial overscan bounding box from an amplifier. 

1191 

1192 This includes the serial/parallel overscan region. 

1193 

1194 Parameters 

1195 ---------- 

1196 amp : `lsst.afw.ampInfo.AmpInfoRecord` 

1197 Amplifier to operate on. 

1198 

1199 Returns 

1200 ------- 

1201 bbox : `lsst.geom.Box2I` 

1202 """ 

1203 # This only works for untrimmed data. 

1204 bbox = amp.getRawDataBBox() 

1205 

1206 parallelOverscanBBox = amp.getRawParallelOverscanBBox() 

1207 grownImageBBox = bbox.expandedTo(parallelOverscanBBox) 

1208 

1209 serialOverscanBBox = amp.getRawSerialOverscanBBox() 

1210 # Extend the serial overscan bbox to include corners 

1211 serialOverscanBBox = geom.Box2I( 

1212 geom.Point2I(serialOverscanBBox.getMinX(), 

1213 grownImageBBox.getMinY()), 

1214 geom.Extent2I(serialOverscanBBox.getWidth(), 

1215 grownImageBBox.getHeight())) 

1216 

1217 return serialOverscanBBox 

1218 

1219 

1220class RawMockLSST(IsrMockLSST): 

1221 """Generate a raw exposure suitable for ISR. 

1222 """ 

1223 def __init__(self, **kwargs): 

1224 super().__init__(**kwargs) 

1225 self.config.isTrimmed = False 

1226 self.config.doGenerateImage = True 

1227 self.config.doGenerateAmpDict = False 

1228 

1229 # Add astro effects 

1230 self.config.doAddSky = True 

1231 self.config.doAddSource = True 

1232 

1233 # Add optical effects 

1234 self.config.doAddFringe = True 

1235 

1236 # Add instrument effects 

1237 self.config.doAddParallelOverscanRamp = True 

1238 self.config.doAddSerialOverscanRamp = True 

1239 self.config.doAddCrosstalk = True 

1240 self.config.doAddBias = True 

1241 self.config.doAddDark = True 

1242 

1243 self.config.doAddFlat = True 

1244 

1245 

1246class TrimmedRawMockLSST(RawMockLSST): 

1247 """Generate a trimmed raw exposure. 

1248 """ 

1249 def __init__(self, **kwargs): 

1250 super().__init__(**kwargs) 

1251 self.config.isTrimmed = True 

1252 self.config.doAddParallelOverscanRamp = False 

1253 self.config.doAddSerialOverscanRamp = False 

1254 

1255 

1256class CalibratedRawMockLSST(RawMockLSST): 

1257 """Generate a trimmed raw exposure. 

1258 

1259 This represents a "truth" image that can be compared to a 

1260 post-ISR cleaned image. 

1261 """ 

1262 def __init__(self, **kwargs): 

1263 super().__init__(**kwargs) 

1264 self.config.isTrimmed = True 

1265 self.config.doGenerateImage = True 

1266 

1267 self.config.doAddSky = True 

1268 self.config.doAddSource = True 

1269 

1270 self.config.doAddFringe = True 

1271 

1272 self.config.doAddParallelOverscanRamp = False 

1273 self.config.doAddSerialOverscanRamp = False 

1274 self.config.doAddCrosstalk = False 

1275 self.config.doAddBias = False 

1276 self.config.doAdd2DBias = False 

1277 self.config.doAddDark = False 

1278 self.config.doApplyGain = False 

1279 self.config.doAddFlat = False 

1280 self.config.doAddClockInjectedOffset = False 

1281 

1282 self.config.biasLevel = 0.0 

1283 # Assume combined calibrations are made with 16 inputs. 

1284 self.config.readNoise *= 0.25 

1285 

1286 self.config.doRoundAdu = False 

1287 

1288 

1289class ReferenceMockLSST(IsrMockLSST): 

1290 """Parent class for those that make reference calibrations. 

1291 """ 

1292 def __init__(self, **kwargs): 

1293 # If we want the calibration in adu units, we need to apply 

1294 # the gain. Default is electron units, so do not apply the gain. 

1295 doApplyGain = kwargs.pop("adu", False) 

1296 

1297 super().__init__(**kwargs) 

1298 self.config.isTrimmed = True 

1299 self.config.doGenerateImage = True 

1300 

1301 self.config.calibMode = True 

1302 

1303 self.config.doAddSky = False 

1304 self.config.doAddSource = False 

1305 

1306 self.config.doAddFringe = False 

1307 

1308 self.config.doAddParallelOverscanRamp = False 

1309 self.config.doAddSerialOverscanRamp = False 

1310 self.config.doAddCrosstalk = False 

1311 self.config.doAddBias = False 

1312 self.config.doAdd2DBias = False 

1313 self.config.doAddDark = False 

1314 self.config.doApplyGain = doApplyGain 

1315 self.config.doAddFlat = False 

1316 self.config.doAddClockInjectedOffset = False 

1317 

1318 # Reference calibrations are not integerized. 

1319 self.config.doRoundAdu = False 

1320 

1321 

1322# Classes to generate calibration products mocks. 

1323class DarkMockLSST(ReferenceMockLSST): 

1324 """Simulated reference dark calibration. 

1325 """ 

1326 def __init__(self, **kwargs): 

1327 super().__init__(**kwargs) 

1328 self.config.doAddDark = True 

1329 self.config.darkTime = 1.0 

1330 

1331 

1332class BiasMockLSST(ReferenceMockLSST): 

1333 """Simulated combined bias calibration. 

1334 """ 

1335 def __init__(self, **kwargs): 

1336 super().__init__(**kwargs) 

1337 # This is a "2D bias residual" frame which has only 

1338 # the 2D bias in it. 

1339 self.config.doAdd2DBias = True 

1340 

1341 

1342class FlatMockLSST(ReferenceMockLSST): 

1343 """Simulated reference flat calibration. 

1344 """ 

1345 def __init__(self, **kwargs): 

1346 super().__init__(**kwargs) 

1347 self.config.doAddFlat = True 

1348 

1349 

1350class FringeMockLSST(ReferenceMockLSST): 

1351 """Simulated reference fringe calibration. 

1352 """ 

1353 def __init__(self, **kwargs): 

1354 super().__init__(**kwargs) 

1355 self.config.doAddFringe = True 

1356 

1357 

1358class BfKernelMockLSST(IsrMockLSST): 

1359 """Simulated brighter-fatter kernel. 

1360 """ 

1361 def __init__(self, **kwargs): 

1362 super().__init__(**kwargs) 

1363 self.config.doGenerateImage = False 

1364 self.config.doGenerateData = True 

1365 

1366 self.config.doBrighterFatter = True 

1367 self.config.doDefects = False 

1368 self.config.doCrosstalkCoeffs = False 

1369 self.config.doTransmissionCurve = False 

1370 self.config.doLinearizer = False 

1371 

1372 

1373class DeferredChargeMockLSST(IsrMockLSST): 

1374 """Simulated deferred charge calibration. 

1375 """ 

1376 def __init__(self, **kwargs): 

1377 super().__init__(**kwargs) 

1378 self.config.doGenerateImage = False 

1379 self.config.doGenerateData = True 

1380 self.config.doDeferredCharge = True 

1381 self.config.doDefects = False 

1382 self.config.doCrosstalkCoeffs = False 

1383 self.config.doTransmissionCurve = False 

1384 

1385 

1386class DefectMockLSST(IsrMockLSST): 

1387 """Simulated defect list. 

1388 """ 

1389 def __init__(self, **kwargs): 

1390 super().__init__(**kwargs) 

1391 self.config.doGenerateImage = False 

1392 self.config.doGenerateData = True 

1393 

1394 self.config.doBrighterFatter = False 

1395 self.config.doDefects = True 

1396 self.config.doCrosstalkCoeffs = False 

1397 self.config.doTransmissionCurve = False 

1398 self.config.doLinearizer = False 

1399 

1400 

1401class CrosstalkCoeffMockLSST(IsrMockLSST): 

1402 """Simulated crosstalk coefficient matrix. 

1403 """ 

1404 def __init__(self, **kwargs): 

1405 super().__init__(**kwargs) 

1406 self.config.doGenerateImage = False 

1407 self.config.doGenerateData = True 

1408 

1409 self.config.doBrighterFatter = False 

1410 self.config.doDefects = False 

1411 self.config.doCrosstalkCoeffs = True 

1412 self.config.doTransmissionCurve = False 

1413 self.config.doLinearizer = False 

1414 

1415 

1416class LinearizerMockLSST(IsrMockLSST): 

1417 """Simulated linearizer. 

1418 """ 

1419 def __init__(self, **kwargs): 

1420 super().__init__(**kwargs) 

1421 self.config.doGenerateImage = False 

1422 self.config.doGenerateData = True 

1423 

1424 self.config.doBrighterFatter = False 

1425 self.config.doDefects = False 

1426 self.config.doCrosstalkCoeffs = False 

1427 self.config.doTransmissionCurve = False 

1428 self.config.doLinearizer = True 

1429 

1430 

1431class TransmissionMockLSST(IsrMockLSST): 

1432 """Simulated transmission curve. 

1433 """ 

1434 def __init__(self, **kwargs): 

1435 super().__init__(**kwargs) 

1436 self.config.doGenerateImage = False 

1437 self.config.doGenerateData = True 

1438 

1439 self.config.doBrighterFatter = False 

1440 self.config.doDefects = False 

1441 self.config.doCrosstalkCoeffs = False 

1442 self.config.doTransmissionCurve = True 

1443 self.config.doLinearizer = False