lsst.ip.isr g2d4321ec69+d787c9b750
ptcDataset.py
Go to the documentation of this file.
2# LSST Data Management System
3# Copyright 2008-2017 AURA/LSST.
4#
5# This product includes software developed by the
6# LSST Project (http://www.lsst.org/).
7#
8# This program is free software: you can redistribute it and/or modify
9# it under the terms of the GNU General Public License as published by
10# the Free Software Foundation, either version 3 of the License, or
11# (at your option) any later version.
12#
13# This program is distributed in the hope that it will be useful,
14# but WITHOUT ANY WARRANTY; without even the implied warranty of
15# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16# GNU General Public License for more details.
17#
18# You should have received a copy of the LSST License Statement and
19# the GNU General Public License along with this program. If not,
20# see <https://www.lsstcorp.org/LegalNotices/>.
21#
22"""
23Define dataset class for MeasurePhotonTransferCurve task
24"""
25import numpy as np
26from astropy.table import Table
27
28from lsst.ip.isr import IsrCalib
29
30__all__ = ['PhotonTransferCurveDataset']
31
32
34 """A simple class to hold the output data from the PTC task.
35 The dataset is made up of a dictionary for each item, keyed by the
36 amplifiers' names, which much be supplied at construction time.
37 New items cannot be added to the class to save accidentally saving to the
38 wrong property, and the class can be frozen if desired.
39 inputExpIdPairs records the exposures used to produce the data.
40 When fitPtc() or fitCovariancesAstier() is run, a mask is built up, which
41 is by definition always the same length as inputExpIdPairs, rawExpTimes,
42 rawMeans and rawVars, and is a list of bools, which are incrementally set
43 to False as points are discarded from the fits.
44 PTC fit parameters for polynomials are stored in a list in ascending order
45 of polynomial term, i.e. par[0]*x^0 + par[1]*x + par[2]*x^2 etc
46 with the length of the list corresponding to the order of the polynomial
47 plus one.
48
49 Parameters
50 ----------
51 ampNames : `list`
52 List with the names of the amplifiers of the detector at hand.
53
54 ptcFitType : `str`
55 Type of model fitted to the PTC: "POLYNOMIAL", "EXPAPPROXIMATION",
56 or "FULLCOVARIANCE".
57
58 covMatrixSide : `int`
59 Maximum lag of covariances (size of square covariance matrices).
60
61 kwargs : `dict`, optional
62 Other keyword arguments to pass to the parent init.
63
64 Notes
65 -----
66 The stored attributes are:
67 badAmps : `list`
68 List with bad amplifiers names.
69 inputExpIdPairs : `dict`, [`str`, `list`]
70 Dictionary keyed by amp names containing the input exposures IDs.
71 expIdMask : `dict`, [`str`, `list`]
72 Dictionary keyed by amp names containing the mask produced after
73 outlier rejection. The mask produced by the "FULLCOVARIANCE"
74 option may differ from the one produced in the other two PTC
75 fit types.
76 rawExpTimes : `dict`, [`str`, `list`]
77 Dictionary keyed by amp names containing the unmasked exposure times.
78 rawMeans : `dict`, [`str`, `list`]
79 Dictionary keyed by amp namescontaining the unmasked average of the
80 means of the exposures in each flat pair.
81 rawVars : `dict`, [`str`, `list`]
82 Dictionary keyed by amp names containing the variance of the
83 difference image of the exposures in each flat pair.
84 gain : `dict`, [`str`, `list`]
85 Dictionary keyed by amp names containing the fitted gains.
86 gainErr : `dict`, [`str`, `list`]
87 Dictionary keyed by amp names containing the errors on the
88 fitted gains.
89 noise : `dict`, [`str`, `list`]
90 Dictionary keyed by amp names containing the fitted noise.
91 noiseErr : `dict`, [`str`, `list`]
92 Dictionary keyed by amp names containing the errors on the fitted
93 noise.
94 ptcFitPars : `dict`, [`str`, `list`]
95 Dictionary keyed by amp names containing the fitted parameters of the
96 PTC model for ptcFitTye in ["POLYNOMIAL", "EXPAPPROXIMATION"].
97 ptcFitParsError : `dict`, [`str`, `list`]
98 Dictionary keyed by amp names containing the errors on the fitted
99 parameters of the PTC model for ptcFitTye in
100 ["POLYNOMIAL", "EXPAPPROXIMATION"].
101 ptcFitChiSq : `dict`, [`str`, `list`]
102 Dictionary keyed by amp names containing the reduced chi squared
103 of the fit for ptcFitTye in ["POLYNOMIAL", "EXPAPPROXIMATION"].
104 ptcTurnoff : `float`
105 Flux value (in ADU) where the variance of the PTC curve starts
106 decreasing consistently.
107 covariances : `dict`, [`str`, `list`]
108 Dictionary keyed by amp names containing a list of measured
109 covariances per mean flux.
110 covariancesModel : `dict`, [`str`, `list`]
111 Dictionary keyed by amp names containinging covariances model
112 (Eq. 20 of Astier+19) per mean flux.
113 covariancesSqrtWeights : `dict`, [`str`, `list`]
114 Dictionary keyed by amp names containinging sqrt. of covariances
115 weights.
116 aMatrix : `dict`, [`str`, `list`]
117 Dictionary keyed by amp names containing the "a" parameters from
118 the model in Eq. 20 of Astier+19.
119 bMatrix : `dict`, [`str`, `list`]
120 Dictionary keyed by amp names containing the "b" parameters from
121 the model in Eq. 20 of Astier+19.
122 covariancesModelNoB : `dict`, [`str`, `list`]
123 Dictionary keyed by amp names containing covariances model
124 (with 'b'=0 in Eq. 20 of Astier+19)
125 per mean flux.
126 aMatrixNoB : `dict`, [`str`, `list`]
127 Dictionary keyed by amp names containing the "a" parameters from the
128 model in Eq. 20 of Astier+19
129 (and 'b' = 0).
130 finalVars : `dict`, [`str`, `list`]
131 Dictionary keyed by amp names containing the masked variance of the
132 difference image of each flat
133 pair. If needed, each array will be right-padded with
134 np.nan to match the length of rawExpTimes.
135 finalModelVars : `dict`, [`str`, `list`]
136 Dictionary keyed by amp names containing the masked modeled
137 variance of the difference image of each flat pair. If needed, each
138 array will be right-padded with np.nan to match the length of
139 rawExpTimes.
140 finalMeans : `dict`, [`str`, `list`]
141 Dictionary keyed by amp names containing the masked average of the
142 means of the exposures in each flat pair. If needed, each array
143 will be right-padded with np.nan to match the length of
144 rawExpTimes.
145 photoCharge : `dict`, [`str`, `list`]
146 Dictionary keyed by amp names containing the integrated photocharge
147 for linearity calibration.
148
149 Returns
150 -------
151 `lsst.cp.pipe.ptc.PhotonTransferCurveDataset`
152 Output dataset from MeasurePhotonTransferCurveTask.
153 """
154
155 _OBSTYPE = 'PTC'
156 _SCHEMA = 'Gen3 Photon Transfer Curve'
157 _VERSION = 1.0
158
159 def __init__(self, ampNames=[], ptcFitType=None, covMatrixSide=1, **kwargs):
160
161 self.ptcFitType = ptcFitType
162 self.ampNames = ampNames
163 self.covMatrixSide = covMatrixSide
164
165 self.badAmps = [np.nan]
166
167 self.inputExpIdPairs = {ampName: [] for ampName in ampNames}
168 self.expIdMask = {ampName: [] for ampName in ampNames}
169 self.rawExpTimes = {ampName: [] for ampName in ampNames}
170 self.rawMeans = {ampName: [] for ampName in ampNames}
171 self.rawVars = {ampName: [] for ampName in ampNames}
172 self.photoCharge = {ampName: [] for ampName in ampNames}
173
174 self.gain = {ampName: np.nan for ampName in ampNames}
175 self.gainErr = {ampName: np.nan for ampName in ampNames}
176 self.noise = {ampName: np.nan for ampName in ampNames}
177 self.noiseErr = {ampName: np.nan for ampName in ampNames}
178
179 self.ptcFitPars = {ampName: [] for ampName in ampNames}
180 self.ptcFitParsError = {ampName: [] for ampName in ampNames}
181 self.ptcFitChiSq = {ampName: np.nan for ampName in ampNames}
182 self.ptcTurnoff = {ampName: np.nan for ampName in ampNames}
183
184 self.covariances = {ampName: [] for ampName in ampNames}
185 self.covariancesModel = {ampName: [] for ampName in ampNames}
186 self.covariancesSqrtWeights = {ampName: [] for ampName in ampNames}
187 self.aMatrix = {ampName: np.nan for ampName in ampNames}
188 self.bMatrix = {ampName: np.nan for ampName in ampNames}
189 self.covariancesModelNoB = {ampName: [] for ampName in ampNames}
190 self.aMatrixNoB = {ampName: np.nan for ampName in ampNames}
191
192 self.finalVars = {ampName: [] for ampName in ampNames}
193 self.finalModelVars = {ampName: [] for ampName in ampNames}
194 self.finalMeans = {ampName: [] for ampName in ampNames}
195
196 super().__init__(**kwargs)
197 self.requiredAttributesrequiredAttributesrequiredAttributes.update(['badAmps', 'inputExpIdPairs', 'expIdMask', 'rawExpTimes',
198 'rawMeans', 'rawVars', 'gain', 'gainErr', 'noise', 'noiseErr',
199 'ptcFitPars', 'ptcFitParsError', 'ptcFitChiSq', 'ptcTurnoff',
200 'aMatrixNoB', 'covariances', 'covariancesModel',
201 'covariancesSqrtWeights', 'covariancesModelNoB',
202 'aMatrix', 'bMatrix', 'finalVars', 'finalModelVars', 'finalMeans',
203 'photoCharge'])
204
205 def setAmpValues(self, ampName, inputExpIdPair=[(np.nan, np.nan)], expIdMask=[np.nan],
206 rawExpTime=[np.nan], rawMean=[np.nan], rawVar=[np.nan], photoCharge=[np.nan],
207 gain=np.nan, gainErr=np.nan, noise=np.nan, noiseErr=np.nan, ptcFitPars=[np.nan],
208 ptcFitParsError=[np.nan], ptcFitChiSq=np.nan, ptcTurnoff=np.nan, covArray=[],
209 covArrayModel=[], covSqrtWeights=[], aMatrix=[], bMatrix=[], covArrayModelNoB=[],
210 aMatrixNoB=[], finalVar=[np.nan], finalModelVar=[np.nan], finalMean=[np.nan]):
211 """Function to initialize an amp of a PhotonTransferCurveDataset.
212
213 Notes
214 -----
215 The parameters are all documented in `init`.
216 """
217 nanMatrix = np.full((self.covMatrixSide, self.covMatrixSide), np.nan)
218 if len(covArray) == 0:
219 covArray = [nanMatrix]
220 if len(covArrayModel) == 0:
221 covArrayModel = [nanMatrix]
222 if len(covSqrtWeights) == 0:
223 covSqrtWeights = [nanMatrix]
224 if len(covArrayModelNoB) == 0:
225 covArrayModelNoB = [nanMatrix]
226 if len(aMatrix) == 0:
227 aMatrix = nanMatrix
228 if len(bMatrix) == 0:
229 bMatrix = nanMatrix
230 if len(aMatrixNoB) == 0:
231 aMatrixNoB = nanMatrix
232
233 self.inputExpIdPairs[ampName] = inputExpIdPair
234 self.expIdMask[ampName] = expIdMask
235 self.rawExpTimes[ampName] = rawExpTime
236 self.rawMeans[ampName] = rawMean
237 self.rawVars[ampName] = rawVar
238 self.photoCharge[ampName] = photoCharge
239 self.gain[ampName] = gain
240 self.gainErr[ampName] = gainErr
241 self.noise[ampName] = noise
242 self.noiseErr[ampName] = noiseErr
243 self.ptcFitPars[ampName] = ptcFitPars
244 self.ptcFitParsError[ampName] = ptcFitParsError
245 self.ptcFitChiSq[ampName] = ptcFitChiSq
246 self.ptcTurnoff[ampName] = ptcTurnoff
247 self.covariances[ampName] = covArray
248 self.covariancesSqrtWeights[ampName] = covSqrtWeights
249 self.covariancesModel[ampName] = covArrayModel
250 self.covariancesModelNoB[ampName] = covArrayModelNoB
251 self.aMatrix[ampName] = aMatrix
252 self.bMatrix[ampName] = bMatrix
253 self.aMatrixNoB[ampName] = aMatrixNoB
254 self.ptcFitPars[ampName] = ptcFitPars
255 self.ptcFitParsError[ampName] = ptcFitParsError
256 self.ptcFitChiSq[ampName] = ptcFitChiSq
257 self.finalVars[ampName] = finalVar
258 self.finalModelVars[ampName] = finalModelVar
259 self.finalMeans[ampName] = finalMean
260
261 def updateMetadata(self, setDate=False, **kwargs):
262 """Update calibration metadata.
263 This calls the base class's method after ensuring the required
264 calibration keywords will be saved.
265 Parameters
266 ----------
267 setDate : `bool`, optional
268 Update the CALIBDATE fields in the metadata to the current
269 time. Defaults to False.
270 kwargs :
271 Other keyword parameters to set in the metadata.
272 """
273 kwargs['PTC_FIT_TYPE'] = self.ptcFitType
274
275 super().updateMetadata(setDate=setDate, **kwargs)
276
277 @classmethod
278 def fromDict(cls, dictionary):
279 """Construct a calibration from a dictionary of properties.
280 Must be implemented by the specific calibration subclasses.
281 Parameters
282 ----------
283 dictionary : `dict`
284 Dictionary of properties.
285 Returns
286 -------
287 calib : `lsst.ip.isr.CalibType`
288 Constructed calibration.
289 Raises
290 ------
291 RuntimeError :
292 Raised if the supplied dictionary is for a different
293 calibration.
294 """
295 calib = cls()
296 if calib._OBSTYPE != dictionary['metadata']['OBSTYPE']:
297 raise RuntimeError(f"Incorrect Photon Transfer Curve dataset supplied. "
298 f"Expected {calib._OBSTYPE}, found {dictionary['metadata']['OBSTYPE']}")
299 calib.setMetadata(dictionary['metadata'])
300 calib.ptcFitType = dictionary['ptcFitType']
301 calib.covMatrixSide = dictionary['covMatrixSide']
302 calib.badAmps = np.array(dictionary['badAmps'], 'str').tolist()
303 # The cov matrices are square
304 covMatrixSide = calib.covMatrixSide
305 # Number of final signal levels
306 covDimensionsProduct = len(np.array(list(dictionary['covariances'].values())[0]).ravel())
307 nSignalPoints = int(covDimensionsProduct/(covMatrixSide*covMatrixSide))
308
309 for ampName in dictionary['ampNames']:
310 calib.ampNames.append(ampName)
311 calib.inputExpIdPairs[ampName] = np.array(dictionary['inputExpIdPairs'][ampName]).tolist()
312 calib.expIdMask[ampName] = np.array(dictionary['expIdMask'][ampName]).tolist()
313 calib.rawExpTimes[ampName] = np.array(dictionary['rawExpTimes'][ampName]).tolist()
314 calib.rawMeans[ampName] = np.array(dictionary['rawMeans'][ampName]).tolist()
315 calib.rawVars[ampName] = np.array(dictionary['rawVars'][ampName]).tolist()
316 calib.gain[ampName] = np.array(dictionary['gain'][ampName]).tolist()
317 calib.gainErr[ampName] = np.array(dictionary['gainErr'][ampName]).tolist()
318 calib.noise[ampName] = np.array(dictionary['noise'][ampName]).tolist()
319 calib.noiseErr[ampName] = np.array(dictionary['noiseErr'][ampName]).tolist()
320 calib.ptcFitPars[ampName] = np.array(dictionary['ptcFitPars'][ampName]).tolist()
321 calib.ptcFitParsError[ampName] = np.array(dictionary['ptcFitParsError'][ampName]).tolist()
322 calib.ptcFitChiSq[ampName] = np.array(dictionary['ptcFitChiSq'][ampName]).tolist()
323 calib.ptcTurnoff[ampName] = np.array(dictionary['ptcTurnoff'][ampName]).tolist()
324 calib.covariances[ampName] = np.array(dictionary['covariances'][ampName]).reshape(
325 (nSignalPoints, covMatrixSide, covMatrixSide)).tolist()
326 calib.covariancesModel[ampName] = np.array(
327 dictionary['covariancesModel'][ampName]).reshape(
328 (nSignalPoints, covMatrixSide, covMatrixSide)).tolist()
329 calib.covariancesSqrtWeights[ampName] = np.array(
330 dictionary['covariancesSqrtWeights'][ampName]).reshape(
331 (nSignalPoints, covMatrixSide, covMatrixSide)).tolist()
332 calib.aMatrix[ampName] = np.array(dictionary['aMatrix'][ampName]).reshape(
333 (covMatrixSide, covMatrixSide)).tolist()
334 calib.bMatrix[ampName] = np.array(dictionary['bMatrix'][ampName]).reshape(
335 (covMatrixSide, covMatrixSide)).tolist()
336 calib.covariancesModelNoB[ampName] = np.array(
337 dictionary['covariancesModelNoB'][ampName]).reshape(
338 (nSignalPoints, covMatrixSide, covMatrixSide)).tolist()
339 calib.aMatrixNoB[ampName] = np.array(
340 dictionary['aMatrixNoB'][ampName]).reshape((covMatrixSide, covMatrixSide)).tolist()
341 calib.finalVars[ampName] = np.array(dictionary['finalVars'][ampName]).tolist()
342 calib.finalModelVars[ampName] = np.array(dictionary['finalModelVars'][ampName]).tolist()
343 calib.finalMeans[ampName] = np.array(dictionary['finalMeans'][ampName]).tolist()
344 calib.photoCharge[ampName] = np.array(dictionary['photoCharge'][ampName]).tolist()
345 calib.updateMetadata()
346 return calib
347
348 def toDict(self):
349 """Return a dictionary containing the calibration properties.
350 The dictionary should be able to be round-tripped through
351 `fromDict`.
352 Returns
353 -------
354 dictionary : `dict`
355 Dictionary of properties.
356 """
358
359 outDict = dict()
360 metadata = self.getMetadata()
361 outDict['metadata'] = metadata
362
363 outDict['ptcFitType'] = self.ptcFitType
364 outDict['covMatrixSide'] = self.covMatrixSide
365 outDict['ampNames'] = self.ampNames
366 outDict['badAmps'] = self.badAmps
367 outDict['inputExpIdPairs'] = self.inputExpIdPairs
368 outDict['expIdMask'] = self.expIdMask
369 outDict['rawExpTimes'] = self.rawExpTimes
370 outDict['rawMeans'] = self.rawMeans
371 outDict['rawVars'] = self.rawVars
372 outDict['gain'] = self.gain
373 outDict['gainErr'] = self.gainErr
374 outDict['noise'] = self.noise
375 outDict['noiseErr'] = self.noiseErr
376 outDict['ptcFitPars'] = self.ptcFitPars
377 outDict['ptcFitParsError'] = self.ptcFitParsError
378 outDict['ptcFitChiSq'] = self.ptcFitChiSq
379 outDict['ptcTurnoff'] = self.ptcTurnoff
380 outDict['covariances'] = self.covariances
381 outDict['covariancesModel'] = self.covariancesModel
382 outDict['covariancesSqrtWeights'] = self.covariancesSqrtWeights
383 outDict['aMatrix'] = self.aMatrix
384 outDict['bMatrix'] = self.bMatrix
385 outDict['covariancesModelNoB'] = self.covariancesModelNoB
386 outDict['aMatrixNoB'] = self.aMatrixNoB
387 outDict['finalVars'] = self.finalVars
388 outDict['finalModelVars'] = self.finalModelVars
389 outDict['finalMeans'] = self.finalMeans
390 outDict['photoCharge'] = self.photoCharge
391
392 return outDict
393
394 @classmethod
395 def fromTable(cls, tableList):
396 """Construct calibration from a list of tables.
397 This method uses the `fromDict` method to create the
398 calibration, after constructing an appropriate dictionary from
399 the input tables.
400 Parameters
401 ----------
402 tableList : `list` [`lsst.afw.table.Table`]
403 List of tables to use to construct the datasetPtc.
404 Returns
405 -------
406 calib : `lsst.cp.pipe.`
407 The calibration defined in the tables.
408 """
409 ptcTable = tableList[0]
410
411 metadata = ptcTable.meta
412 inDict = dict()
413 inDict['metadata'] = metadata
414 inDict['ampNames'] = []
415 inDict['ptcFitType'] = []
416 inDict['covMatrixSide'] = []
417 inDict['inputExpIdPairs'] = dict()
418 inDict['expIdMask'] = dict()
419 inDict['rawExpTimes'] = dict()
420 inDict['rawMeans'] = dict()
421 inDict['rawVars'] = dict()
422 inDict['gain'] = dict()
423 inDict['gainErr'] = dict()
424 inDict['noise'] = dict()
425 inDict['noiseErr'] = dict()
426 inDict['ptcFitPars'] = dict()
427 inDict['ptcFitParsError'] = dict()
428 inDict['ptcFitChiSq'] = dict()
429 inDict['ptcTurnoff'] = dict()
430 inDict['covariances'] = dict()
431 inDict['covariancesModel'] = dict()
432 inDict['covariancesSqrtWeights'] = dict()
433 inDict['aMatrix'] = dict()
434 inDict['bMatrix'] = dict()
435 inDict['covariancesModelNoB'] = dict()
436 inDict['aMatrixNoB'] = dict()
437 inDict['finalVars'] = dict()
438 inDict['finalModelVars'] = dict()
439 inDict['finalMeans'] = dict()
440 inDict['badAmps'] = []
441 inDict['photoCharge'] = dict()
442
443 for record in ptcTable:
444 ampName = record['AMPLIFIER_NAME']
445
446 inDict['ptcFitType'] = record['PTC_FIT_TYPE']
447 inDict['covMatrixSide'] = record['COV_MATRIX_SIDE']
448 inDict['ampNames'].append(ampName)
449 inDict['inputExpIdPairs'][ampName] = record['INPUT_EXP_ID_PAIRS']
450 inDict['expIdMask'][ampName] = record['EXP_ID_MASK']
451 inDict['rawExpTimes'][ampName] = record['RAW_EXP_TIMES']
452 inDict['rawMeans'][ampName] = record['RAW_MEANS']
453 inDict['rawVars'][ampName] = record['RAW_VARS']
454 inDict['gain'][ampName] = record['GAIN']
455 inDict['gainErr'][ampName] = record['GAIN_ERR']
456 inDict['noise'][ampName] = record['NOISE']
457 inDict['noiseErr'][ampName] = record['NOISE_ERR']
458 inDict['ptcFitPars'][ampName] = record['PTC_FIT_PARS']
459 inDict['ptcFitParsError'][ampName] = record['PTC_FIT_PARS_ERROR']
460 inDict['ptcFitChiSq'][ampName] = record['PTC_FIT_CHI_SQ']
461 inDict['ptcTurnoff'][ampName] = record['PTC_TURNOFF']
462 inDict['covariances'][ampName] = record['COVARIANCES']
463 inDict['covariancesModel'][ampName] = record['COVARIANCES_MODEL']
464 inDict['covariancesSqrtWeights'][ampName] = record['COVARIANCES_SQRT_WEIGHTS']
465 inDict['aMatrix'][ampName] = record['A_MATRIX']
466 inDict['bMatrix'][ampName] = record['B_MATRIX']
467 inDict['covariancesModelNoB'][ampName] = record['COVARIANCES_MODEL_NO_B']
468 inDict['aMatrixNoB'][ampName] = record['A_MATRIX_NO_B']
469 inDict['finalVars'][ampName] = record['FINAL_VARS']
470 inDict['finalModelVars'][ampName] = record['FINAL_MODEL_VARS']
471 inDict['finalMeans'][ampName] = record['FINAL_MEANS']
472 inDict['badAmps'] = record['BAD_AMPS']
473 inDict['photoCharge'][ampName] = record['PHOTO_CHARGE']
474 return cls().fromDict(inDict)
475
476 def toTable(self):
477 """Construct a list of tables containing the information in this
478 calibration.
479
480 The list of tables should create an identical calibration
481 after being passed to this class's fromTable method.
482 Returns
483 -------
484 tableList : `list` [`astropy.table.Table`]
485 List of tables containing the linearity calibration
486 information.
487 """
488 tableList = []
490 nPoints = []
491 for i, ampName in enumerate(self.ampNames):
492 nPoints.append(len(list(self.covariances.values())[i]))
493 nSignalPoints = max(nPoints)
494 nPadPoints = {}
495 for i, ampName in enumerate(self.ampNames):
496 nPadPoints[ampName] = nSignalPoints - len(list(self.covariances.values())[i])
497 covMatrixSide = self.covMatrixSide
498
499 catalog = Table([{'AMPLIFIER_NAME': ampName,
500 'PTC_FIT_TYPE': self.ptcFitType,
501 'COV_MATRIX_SIDE': self.covMatrixSide,
502 'INPUT_EXP_ID_PAIRS': self.inputExpIdPairs[ampName]
503 if len(self.expIdMask[ampName]) else np.nan,
504 'EXP_ID_MASK': self.expIdMask[ampName]
505 if len(self.expIdMask[ampName]) else np.nan,
506 'RAW_EXP_TIMES': np.array(self.rawExpTimes[ampName]).tolist()
507 if len(self.rawExpTimes[ampName]) else np.nan,
508 'RAW_MEANS': np.array(self.rawMeans[ampName]).tolist()
509 if len(self.rawMeans[ampName]) else np.nan,
510 'RAW_VARS': np.array(self.rawVars[ampName]).tolist()
511 if len(self.rawVars[ampName]) else np.nan,
512 'GAIN': self.gain[ampName],
513 'GAIN_ERR': self.gainErr[ampName],
514 'NOISE': self.noise[ampName],
515 'NOISE_ERR': self.noiseErr[ampName],
516 'PTC_FIT_PARS': np.array(self.ptcFitPars[ampName]).tolist(),
517 'PTC_FIT_PARS_ERROR': np.array(self.ptcFitParsError[ampName]).tolist(),
518 'PTC_FIT_CHI_SQ': self.ptcFitChiSq[ampName],
519 'PTC_TURNOFF': self.ptcTurnoff[ampName],
520 'COVARIANCES': np.pad(np.array(self.covariances[ampName]),
521 ((0, nPadPoints[ampName]), (0, 0), (0, 0)),
522 'constant', constant_values=np.nan).reshape(
523 nSignalPoints*covMatrixSide**2).tolist(),
524 'COVARIANCES_MODEL': np.pad(np.array(self.covariancesModel[ampName]),
525 ((0, nPadPoints[ampName]), (0, 0), (0, 0)),
526 'constant', constant_values=np.nan).reshape(
527 nSignalPoints*covMatrixSide**2).tolist(),
528 'COVARIANCES_SQRT_WEIGHTS': np.pad(np.array(self.covariancesSqrtWeights[ampName]),
529 ((0, nPadPoints[ampName]), (0, 0), (0, 0)),
530 'constant', constant_values=0.0).reshape(
531 nSignalPoints*covMatrixSide**2).tolist(),
532 'A_MATRIX': np.array(self.aMatrix[ampName]).reshape(covMatrixSide**2).tolist(),
533 'B_MATRIX': np.array(self.bMatrix[ampName]).reshape(covMatrixSide**2).tolist(),
534 'COVARIANCES_MODEL_NO_B':
535 np.pad(np.array(self.covariancesModelNoB[ampName]),
536 ((0, nPadPoints[ampName]), (0, 0), (0, 0)),
537 'constant', constant_values=np.nan).reshape(
538 nSignalPoints*covMatrixSide**2).tolist(),
539 'A_MATRIX_NO_B': np.array(self.aMatrixNoB[ampName]).reshape(
540 covMatrixSide**2).tolist(),
541 'FINAL_VARS': np.pad(np.array(self.finalVars[ampName]), (0, nPadPoints[ampName]),
542 'constant', constant_values=np.nan).tolist(),
543 'FINAL_MODEL_VARS': np.pad(np.array(self.finalModelVars[ampName]),
544 (0, nPadPoints[ampName]),
545 'constant', constant_values=np.nan).tolist(),
546 'FINAL_MEANS': np.pad(np.array(self.finalMeans[ampName]),
547 (0, nPadPoints[ampName]),
548 'constant', constant_values=np.nan).tolist(),
549 'BAD_AMPS': np.array(self.badAmps).tolist() if len(self.badAmps) else np.nan,
550 'PHOTO_CHARGE': np.array(self.photoCharge[ampName]).tolist(),
551 } for ampName in self.ampNames])
552 inMeta = self.getMetadata().toDict()
553 outMeta = {k: v for k, v in inMeta.items() if v is not None}
554 outMeta.update({k: "" for k, v in inMeta.items() if v is None})
555 catalog.meta = outMeta
556 tableList.append(catalog)
557
558 return(tableList)
559
560 def getExpIdsUsed(self, ampName):
561 """Get the exposures used, i.e. not discarded, for a given amp.
562 If no mask has been created yet, all exposures are returned.
563 """
564 if len(self.expIdMask[ampName]) == 0:
565 return self.inputExpIdPairs[ampName]
566
567 # if the mask exists it had better be the same length as the expIdPairs
568 assert len(self.expIdMask[ampName]) == len(self.inputExpIdPairs[ampName])
569
570 pairs = self.inputExpIdPairs[ampName]
571 mask = self.expIdMask[ampName]
572 # cast to bool required because numpy
573 return [(exp1, exp2) for ((exp1, exp2), m) in zip(pairs, mask) if bool(m) is True]
574
575 def getGoodAmps(self):
576 return [amp for amp in self.ampNames if amp not in self.badAmps]
def requiredAttributes(self, value)
Definition: calibType.py:142
def updateMetadata(self, camera=None, detector=None, filterName=None, setCalibId=False, setCalibInfo=False, setDate=False, **kwargs)
Definition: calibType.py:181
def __init__(self, ampNames=[], ptcFitType=None, covMatrixSide=1, **kwargs)
Definition: ptcDataset.py:159
def updateMetadata(self, setDate=False, **kwargs)
Definition: ptcDataset.py:261
def setAmpValues(self, ampName, inputExpIdPair=[(np.nan, np.nan)], expIdMask=[np.nan], rawExpTime=[np.nan], rawMean=[np.nan], rawVar=[np.nan], photoCharge=[np.nan], gain=np.nan, gainErr=np.nan, noise=np.nan, noiseErr=np.nan, ptcFitPars=[np.nan], ptcFitParsError=[np.nan], ptcFitChiSq=np.nan, ptcTurnoff=np.nan, covArray=[], covArrayModel=[], covSqrtWeights=[], aMatrix=[], bMatrix=[], covArrayModelNoB=[], aMatrixNoB=[], finalVar=[np.nan], finalModelVar=[np.nan], finalMean=[np.nan])
Definition: ptcDataset.py:210