336 inputExpIdPair=(-1, -1),
340 rowMeanVariance=np.nan,
353 Set the amp values for a partial PTC Dataset (from cpExtractPtcTask).
358 Name of the amp to set the values.
359 inputExpIdPair : `tuple` [`int`]
360 Exposure IDs of input pair.
361 rawExpTime : `float`, optional
362 Exposure time for this exposure pair.
363 rawMean : `float`, optional
364 Average of the means of the exposures in this pair.
365 rawVar : `float`, optional
366 Variance of the difference of the exposures in this pair.
367 rowMeanVariance : `float`, optional
368 Variance of the means of the rows in the difference image
369 of the exposures in this pair.
370 photoCharge : `float`, optional
371 Integrated photocharge for flat pair for linearity calibration.
372 ampOffset : `float`, optional
373 Amp offset for this amplifier.
374 expIdMask : `bool`, optional
375 Flag setting if this exposure pair should be used (True)
377 covariance : `np.ndarray` or None, optional
378 Measured covariance for this exposure pair.
379 covSqrtWeights : `np.ndarray` or None, optional
380 Measured sqrt of covariance weights in this exposure pair.
381 gain : `float`, optional
382 Estimated gain for this exposure pair.
383 noise : `float`, optional
384 Estimated read noise for this exposure pair.
385 histVar : `float`, optional
386 Variance estimated from fitting a histogram with a Gaussian model.
387 histChi2Dof : `float`, optional
388 Chi-squared per degree of freedom from Gaussian histogram fit.
389 kspValue : `float`, optional
390 KS test p-value from the Gaussian histogram fit.
395 if covariance
is None:
396 covariance = nanMatrix
397 if covSqrtWeights
is None:
398 covSqrtWeights = nanMatrix
402 self.
rawMeans[ampName] = np.array([rawMean])
403 self.
rawVars[ampName] = np.array([rawVar])
406 self.
ampOffsets[ampName] = np.array([ampOffset])
407 self.
expIdMask[ampName] = np.array([expIdMask])
410 self.
gain[ampName] = gain
412 self.
gainList[ampName] = np.array([gain])
413 self.
noise[ampName] = noise
414 self.
noiseList[ampName] = np.array([noise])
415 self.
histVars[ampName] = np.array([histVar])
417 self.
kspValues[ampName] = np.array([kspValue])
421 self.
aMatrix[ampName] = nanMatrixFit
422 self.
bMatrix[ampName] = nanMatrixFit
425 self.
finalVars[ampName] = np.array([np.nan])
468 """Construct a calibration from a dictionary of properties.
469 Must be implemented by the specific calibration subclasses.
474 Dictionary of properties.
478 calib : `lsst.ip.isr.PhotonTransferCurveDataset`
479 Constructed calibration.
484 Raised if the supplied dictionary is for a different
488 if calib._OBSTYPE != dictionary[
'metadata'][
'OBSTYPE']:
489 raise RuntimeError(f
"Incorrect Photon Transfer Curve dataset supplied. "
490 f
"Expected {calib._OBSTYPE}, found {dictionary['metadata']['OBSTYPE']}")
491 calib.setMetadata(dictionary[
'metadata'])
492 calib.ptcFitType = dictionary[
'ptcFitType']
493 calib.covMatrixSide = dictionary[
'covMatrixSide']
494 calib.covMatrixSideFullCovFit = dictionary[
'covMatrixSideFullCovFit']
495 calib.badAmps = np.array(dictionary[
'badAmps'],
'str').tolist()
498 calibVersion = dictionary[
'metadata'][
'PTC_VERSION']
501 covMatrixSide = calib.covMatrixSide
502 covMatrixSideFullCovFit = calib.covMatrixSideFullCovFit
504 covDimensionsProduct = len(np.array(list(dictionary[
'covariances'].values())[0]).ravel())
505 nSignalPoints = int(covDimensionsProduct/(covMatrixSide*covMatrixSide))
507 for ampName
in dictionary[
'ampNames']:
508 calib.ampNames.append(ampName)
509 calib.inputExpIdPairs[ampName] = dictionary[
'inputExpIdPairs'][ampName]
510 calib.expIdMask[ampName] = np.array(dictionary[
'expIdMask'][ampName])
511 calib.rawExpTimes[ampName] = np.array(dictionary[
'rawExpTimes'][ampName], dtype=np.float64)
512 calib.rawMeans[ampName] = np.array(dictionary[
'rawMeans'][ampName], dtype=np.float64)
513 calib.rawVars[ampName] = np.array(dictionary[
'rawVars'][ampName], dtype=np.float64)
514 calib.rowMeanVariance[ampName] = np.array(dictionary[
'rowMeanVariance'][ampName],
516 calib.gain[ampName] = float(dictionary[
'gain'][ampName])
517 calib.gainErr[ampName] = float(dictionary[
'gainErr'][ampName])
518 calib.gainUnadjusted[ampName] = float(dictionary[
'gainUnadjusted'][ampName])
519 calib.gainList[ampName] = np.array(dictionary[
'gainList'][ampName], dtype=np.float64)
520 calib.noiseList[ampName] = np.array(dictionary[
'noiseList'][ampName], dtype=np.float64)
521 calib.noise[ampName] = float(dictionary[
'noise'][ampName])
522 calib.noiseErr[ampName] = float(dictionary[
'noiseErr'][ampName])
523 calib.histVars[ampName] = np.array(dictionary[
'histVars'][ampName], dtype=np.float64)
524 calib.histChi2Dofs[ampName] = np.array(dictionary[
'histChi2Dofs'][ampName], dtype=np.float64)
525 calib.kspValues[ampName] = np.array(dictionary[
'kspValues'][ampName], dtype=np.float64)
526 calib.ptcFitPars[ampName] = np.array(dictionary[
'ptcFitPars'][ampName], dtype=np.float64)
527 calib.ptcFitParsError[ampName] = np.array(dictionary[
'ptcFitParsError'][ampName],
529 calib.ptcFitChiSq[ampName] = float(dictionary[
'ptcFitChiSq'][ampName])
530 calib.ptcTurnoff[ampName] = float(dictionary[
'ptcTurnoff'][ampName])
531 calib.ptcTurnoffSamplingError[ampName] = float(dictionary[
'ptcTurnoffSamplingError'][ampName])
532 if nSignalPoints > 0:
534 calib.covariances[ampName] = np.array(dictionary[
'covariances'][ampName],
535 dtype=np.float64).reshape(
536 (nSignalPoints, covMatrixSide, covMatrixSide))
537 calib.covariancesModel[ampName] = np.array(
538 dictionary[
'covariancesModel'][ampName],
539 dtype=np.float64).reshape(
540 (nSignalPoints, covMatrixSideFullCovFit, covMatrixSideFullCovFit))
541 calib.covariancesSqrtWeights[ampName] = np.array(
542 dictionary[
'covariancesSqrtWeights'][ampName],
543 dtype=np.float64).reshape(
544 (nSignalPoints, covMatrixSide, covMatrixSide))
545 calib.aMatrix[ampName] = np.array(dictionary[
'aMatrix'][ampName],
546 dtype=np.float64).reshape(
547 (covMatrixSideFullCovFit, covMatrixSideFullCovFit))
548 calib.bMatrix[ampName] = np.array(dictionary[
'bMatrix'][ampName],
549 dtype=np.float64).reshape(
550 (covMatrixSideFullCovFit, covMatrixSideFullCovFit))
551 calib.noiseMatrix[ampName] = np.array(
552 dictionary[
'noiseMatrix'][ampName],
553 dtype=np.float64).reshape((covMatrixSideFullCovFit, covMatrixSideFullCovFit))
557 if calibVersion < 2.1:
558 calib._covariancesModelNoB[ampName] = np.array(
559 dictionary[
'covariancesModelNoB'][ampName], dtype=np.float64).reshape(
560 (nSignalPoints, covMatrixSideFullCovFit, covMatrixSideFullCovFit))
561 calib._aMatrixNoB[ampName] = np.array(
562 dictionary[
'aMatrixNoB'][ampName],
563 dtype=np.float64).reshape((covMatrixSideFullCovFit, covMatrixSideFullCovFit))
564 calib._noiseMatrixNoB[ampName] = np.array(
565 dictionary[
'noiseMatrixNoB'][ampName],
566 dtype=np.float64).reshape((covMatrixSideFullCovFit, covMatrixSideFullCovFit))
570 calib.covariances[ampName] = np.array([], dtype=np.float64)
571 calib.covariancesModel[ampName] = np.array([], dtype=np.float64)
572 calib.covariancesSqrtWeights[ampName] = np.array([], dtype=np.float64)
573 calib.aMatrix[ampName] = np.array([], dtype=np.float64)
574 calib.bMatrix[ampName] = np.array([], dtype=np.float64)
575 calib.noiseMatrix[ampName] = np.array([], dtype=np.float64)
577 calib.finalVars[ampName] = np.array(dictionary[
'finalVars'][ampName], dtype=np.float64)
578 calib.finalModelVars[ampName] = np.array(dictionary[
'finalModelVars'][ampName], dtype=np.float64)
579 calib.finalMeans[ampName] = np.array(dictionary[
'finalMeans'][ampName], dtype=np.float64)
580 calib.photoCharges[ampName] = np.array(dictionary[
'photoCharges'][ampName], dtype=np.float64)
581 calib.ampOffsets[ampName] = np.array(dictionary[
'ampOffsets'][ampName], dtype=np.float64)
583 for key, value
in dictionary[
'auxValues'].
items():
584 if isinstance(value[0], numbers.Integral):
585 calib.auxValues[key] = np.atleast_1d(np.asarray(value).astype(np.int64))
586 elif isinstance(value[0], (str, np.str_, np.string_)):
587 calib.auxValues[key] = np.atleast_1d(np.asarray(value))
589 calib.auxValues[key] = np.atleast_1d(np.array(value, dtype=np.float64))
591 calib.updateMetadata()
595 """Return a dictionary containing the calibration properties.
596 The dictionary should be able to be round-tripped through
602 Dictionary of properties.
608 outDict[
'metadata'] = metadata
610 def _dictOfArraysToDictOfLists(dictOfArrays):
612 for key, value
in dictOfArrays.items():
613 dictOfLists[key] = value.ravel().tolist()
621 outDict[
'badAmps'] = self.
badAmps
623 outDict[
'expIdMask'] = _dictOfArraysToDictOfLists(self.
expIdMask)
624 outDict[
'rawExpTimes'] = _dictOfArraysToDictOfLists(self.
rawExpTimes)
625 outDict[
'rawMeans'] = _dictOfArraysToDictOfLists(self.
rawMeans)
626 outDict[
'rawVars'] = _dictOfArraysToDictOfLists(self.
rawVars)
627 outDict[
'rowMeanVariance'] = _dictOfArraysToDictOfLists(self.
rowMeanVariance)
628 outDict[
'gain'] = self.
gain
629 outDict[
'gainErr'] = self.
gainErr
631 outDict[
'gainList'] = _dictOfArraysToDictOfLists(self.
gainList)
632 outDict[
'noiseList'] = _dictOfArraysToDictOfLists(self.
noiseList)
633 outDict[
'noise'] = self.
noise
638 outDict[
'ptcFitPars'] = _dictOfArraysToDictOfLists(self.
ptcFitPars)
639 outDict[
'ptcFitParsError'] = _dictOfArraysToDictOfLists(self.
ptcFitParsError)
643 outDict[
'covariances'] = _dictOfArraysToDictOfLists(self.
covariances)
644 outDict[
'covariancesModel'] = _dictOfArraysToDictOfLists(self.
covariancesModel)
646 outDict[
'aMatrix'] = _dictOfArraysToDictOfLists(self.
aMatrix)
647 outDict[
'bMatrix'] = _dictOfArraysToDictOfLists(self.
bMatrix)
648 outDict[
'noiseMatrix'] = _dictOfArraysToDictOfLists(self.
noiseMatrix)
649 outDict[
'finalVars'] = _dictOfArraysToDictOfLists(self.
finalVars)
650 outDict[
'finalModelVars'] = _dictOfArraysToDictOfLists(self.
finalModelVars)
651 outDict[
'finalMeans'] = _dictOfArraysToDictOfLists(self.
finalMeans)
652 outDict[
'photoCharges'] = _dictOfArraysToDictOfLists(self.
photoCharges)
653 outDict[
'ampOffsets'] = _dictOfArraysToDictOfLists(self.
ampOffsets)
654 outDict[
'auxValues'] = _dictOfArraysToDictOfLists(self.
auxValues)
941 """Append a partial PTC dataset to this dataset.
945 partialPtc : `lsst.ip.isr.PhotonTransferCurveDataset`
946 Partial PTC to append. Should only have one element.
948 if self.ampNames != partialPtc.ampNames:
949 raise ValueError(
"partialPtc has mis-matched amps.")
950 if len(partialPtc.rawMeans[self.ampNames[0]]) != 1
or partialPtc.ptcFitType !=
"PARTIAL":
951 raise ValueError(
"partialPtc does not appear to be the correct format.")
954 initialLength = len(self.expIdMask[self.ampNames[0]])
956 for key, value
in partialPtc.auxValues.items():
957 if key
in self.auxValues:
958 self.auxValues[key] = np.append(self.auxValues[key], value)
959 elif initialLength == 0:
961 self.auxValues[key] = value
963 raise ValueError(f
"partialPtc has mismatched auxValue key {key}.")
965 for ampName
in self.ampNames:
970 self.inputExpIdPairs[ampName].append(partialPtc.inputExpIdPairs[ampName][0])
971 self.expIdMask[ampName] = np.append(self.expIdMask[ampName],
972 partialPtc.expIdMask[ampName][0])
973 self.rawExpTimes[ampName] = np.append(self.rawExpTimes[ampName],
974 partialPtc.rawExpTimes[ampName][0])
975 self.rawMeans[ampName] = np.append(self.rawMeans[ampName],
976 partialPtc.rawMeans[ampName][0])
977 self.rawVars[ampName] = np.append(self.rawVars[ampName],
978 partialPtc.rawVars[ampName][0])
979 self.rowMeanVariance[ampName] = np.append(self.rowMeanVariance[ampName],
980 partialPtc.rowMeanVariance[ampName][0])
981 self.photoCharges[ampName] = np.append(self.photoCharges[ampName],
982 partialPtc.photoCharges[ampName][0])
983 self.ampOffsets[ampName] = np.append(self.ampOffsets[ampName],
984 partialPtc.ampOffsets[ampName][0])
985 self.histVars[ampName] = np.append(self.histVars[ampName],
986 partialPtc.histVars[ampName][0])
987 self.histChi2Dofs[ampName] = np.append(self.histChi2Dofs[ampName],
988 partialPtc.histChi2Dofs[ampName][0])
989 self.kspValues[ampName] = np.append(self.kspValues[ampName],
990 partialPtc.kspValues[ampName][0])
991 self.gainList[ampName] = np.append(self.gainList[ampName],
992 partialPtc.gain[ampName])
993 self.noiseList[ampName] = np.append(self.noiseList[ampName],
994 partialPtc.noise[ampName])
995 self.finalVars[ampName] = np.append(self.finalVars[ampName],
996 partialPtc.finalVars[ampName][0])
997 self.finalModelVars[ampName] = np.append(self.finalModelVars[ampName],
998 partialPtc.finalModelVars[ampName][0])
999 self.finalMeans[ampName] = np.append(self.finalMeans[ampName],
1000 partialPtc.finalMeans[ampName][0])
1002 self.covariances[ampName] = np.append(
1003 self.covariances[ampName].ravel(),
1004 partialPtc.covariances[ampName].ravel()
1007 len(self.rawExpTimes[ampName]),
1012 self.covariancesSqrtWeights[ampName] = np.append(
1013 self.covariancesSqrtWeights[ampName].ravel(),
1014 partialPtc.covariancesSqrtWeights[ampName].ravel()
1017 len(self.rawExpTimes[ampName]),
1022 self.covariancesModel[ampName] = np.append(
1023 self.covariancesModel[ampName].ravel(),
1024 partialPtc.covariancesModel[ampName].ravel()
1027 len(self.rawExpTimes[ampName]),
1183 """Computes the covariance model at specific signal levels.
1187 mu : `numpy.array`, (N,)
1188 List of mean signals in ADU.
1193 Raised if ptcFitType is invalid.
1197 covModel : `numpy.array`, (N, M, M)
1198 Covariances model at mu (in ADU^2).
1202 Computes the covModel for all mu, and it returns
1203 cov[N, M, M], where the variance model is cov[:,0,0].
1204 Both mu and cov are in ADUs and ADUs squared. This
1205 routine evaulates the n-degree polynomial model (defined
1206 by polynomialFitDegree) if self.ptcFitType == POLYNOMIAL,
1207 the approximation in Eq. 16 of Astier+19 (1905.08677)
1208 if self.ptcFitType == EXPAPPROXIMATION, and Eq. 20 of
1209 Astier+19 if self.ptcFitType == FULLCOVARIANCE.
1211 The POLYNOMIAL model and the EXPAPPROXIMATION model
1212 (Eq. 16 of Astier+19) are only approximations for the
1213 variance (cov[0,0]), so the function returns covModel
1214 of shape (N,), representing an array of [C_{00}]
1215 if self.ptcFitType == EXPAPPROXIMATION or
1216 self.ptcFitType == POLYNOMAIL.
1220 covModel = {ampName: np.array([])
for ampName
in ampNames}
1225 for ampName
in ampNames:
1226 c00 = poly.polyval(mu, [*pars[ampName]])
1227 covModel[ampName] = c00
1232 for ampName
in ampNames:
1233 a00, gain, noise = pars[ampName]
1234 f1 = 0.5/(a00*gain*gain)*(np.exp(2*a00*mu*gain)-1)
1235 f2 = noise/(gain*gain)
1237 covModel[ampName] = c00
1239 elif self.
ptcFitType in [
"FULLCOVARIANCE",
"FULLCOVARIANCE_NO_B"]:
1240 for ampName
in ampNames:
1242 gain = self.
gain[ampName]
1243 aMatrix = self.
aMatrix[ampName]
1244 bMatrix = self.
bMatrix[ampName]
1245 cMatrix = aMatrix*bMatrix
1248 sa = (matrixSideFit, matrixSideFit)
1251 aEnlarged = np.zeros((int(sa[0]*1.5)+1, int(sa[1]*1.5)+1))
1252 aEnlarged[0:sa[0], 0:sa[1]] = aMatrix
1256 cEnlarged = np.zeros((int(sa[0]*1.5)+1, int(sa[1]*1.5)+1))
1257 cEnlarged[0:sa[0], 0:sa[1]] = cMatrix
1260 a2 = fftconvolve(aSym, aSym, mode=
'same')
1261 a3 = fftconvolve(a2, aSym, mode=
'same')
1262 ac = fftconvolve(aSym, cSym, mode=
'same')
1263 (xc, yc) = np.unravel_index(np.abs(aSym).argmax(), a2.shape)
1265 a1 = aMatrix[np.newaxis, :, :]
1266 a2 = a2[np.newaxis, xc:xc + matrixSideFit, yc:yc + matrixSideFit]
1267 a3 = a3[np.newaxis, xc:xc + matrixSideFit, yc:yc + matrixSideFit]
1268 ac = ac[np.newaxis, xc:xc + matrixSideFit, yc:yc + matrixSideFit]
1269 c1 = cMatrix[np.newaxis, ::]
1272 bigMu = mu[:, np.newaxis, np.newaxis]*gain
1276 covModel[ampName] = (bigMu/(gain*gain)*(a1*bigMu+2./3.*(bigMu*bigMu)*(a2 + c1)
1277 + (1./3.*a3 + 5./6.*ac)*(bigMu*bigMu*bigMu))
1278 + noiseMatrix[np.newaxis, :, :]/gain**2)
1281 covModel[ampName][:, 0, 0] += mu/gain
1283 raise RuntimeError(
"Cannot compute PTC model for "