Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

DM-42165: eo_pipe/cp_verify parity: bfk #40

Merged
merged 7 commits into from
Jan 31, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
41 changes: 31 additions & 10 deletions python/lsst/cp/verify/repackStats.py
Original file line number Diff line number Diff line change
Expand Up @@ -378,20 +378,41 @@ def repackDetStats(self, detectorStats, detectorDims):
"ptcFitType": stats["PTC_FIT_TYPE"],
"ptcBfeA00": stats["PTC_BFE_A00"],
"ptcRowMeanVariance": stats["PTC_ROW_MEAN_VARIANCE"],
"ptcRowMeanVarianceSlope": stats["PTC_ROW_MEAN_VARIANCE_SLOPE"],
"ptcMaxRawMeans": stats["PTC_MAX_RAW_MEANS"],
"ptcRawMeans": stats["PTC_RAW_MEANS"],
"ptcExpIdmask": stats["PTC_EXP_ID_MASK"],
"ptcCov10": stats["PTC_COV_10"],
"ptcCov10FitSlope": stats["PTC_COV_10_FIT_SLOPE"],
"ptcCov10FitOffset": stats["PTC_COV_10_FIT_OFFSET"],
"ptcCov10FitSuccess": stats["PTC_COV_10_FIT_SUCCESS"],
"ptcCov01": stats["PTC_COV_01"],
"ptcCov01FitSlope": stats["PTC_COV_01_FIT_SLOPE"],
"ptcCov01FitOffset": stats["PTC_COV_01_FIT_OFFSET"],
"ptcCov01FitSuccess": stats["PTC_COV_01_FIT_SUCCESS"],
"ptcCov11": stats["PTC_COV_11"],
"ptcCov11FitSlope": stats["PTC_COV_11_FIT_SLOPE"],
"ptcCov11FitOffset": stats["PTC_COV_11_FIT_OFFSET"],
"ptcCov11FitSuccess": stats["PTC_COV_11_FIT_SUCCESS"],
"ptcCov20": stats["PTC_COV_20"],
"ptcCov20FitSlope": stats["PTC_COV_20_FIT_SLOPE"],
"ptcCov20FitOffset": stats["PTC_COV_20_FIT_OFFSET"],
"ptcCov20FitSuccess": stats["PTC_COV_20_FIT_SUCCESS"],
"ptcCov02": stats["PTC_COV_02"],
"ptcCov02FitSlope": stats["PTC_COV_02_FIT_SLOPE"],
"ptcCov02FitOffset": stats["PTC_COV_02_FIT_OFFSET"],
"ptcCov02FitSuccess": stats["PTC_COV_02_FIT_SUCCESS"],
}
# Get catalog stats
# Get detector stats
# Get metadata stats
# Get verify stats; no need to loop here.
stats = detStats["VERIFY"]
row["detector"] = {
"instrument": instrument,
"detector": detector,
"ptcVerifyGain": stats["PTC_GAIN"],
"ptcVerifyNoise": stats["PTC_NOISE"],
"ptcVerifyTurnoff": stats["PTC_TURNOFF"],
"ptcVerifyBfeA00": stats["PTC_BFE_A00"],
}
# Get verify stats
for ampName, stats in detStats["VERIFY"]["AMP"].items():
row[ampName]["ptcVerifyGain"] = stats["PTC_GAIN"]
row[ampName]["ptcVerifyNoise"] = stats["PTC_NOISE"]
row[ampName]["ptcVerifyTurnoff"] = stats["PTC_TURNOFF"]
row[ampName]["ptcVerifyBfeA00"] = stats["PTC_BFE_A00"]

# Get isr stats

# Append to output
Expand Down
64 changes: 42 additions & 22 deletions python/lsst/cp/verify/verifyPtc.py
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,7 @@
# along with this program. If not, see <http://www.gnu.org/licenses/>.
import numpy as np
import lsst.pex.config as pexConfig
from scipy.optimize import least_squares

from .verifyCalib import CpVerifyCalibConfig, CpVerifyCalibTask, CpVerifyCalibConnections

Expand Down Expand Up @@ -78,6 +79,18 @@ def setDefaults(self):
)


def linearModel(x, m, b):
"""A linear model.
"""
return m*x + b


def modelResidual(p, x, y):
"""Model residual for fit below.
"""
return y - linearModel(x, *p)


class CpVerifyPtcTask(CpVerifyCalibTask):
"""PTC verification sub-class, implementing the verify method.
"""
Expand Down Expand Up @@ -132,6 +145,33 @@ def amplifierStatistics(self, inputCalib, camera=None):
outputStatistics[ampName]['PTC_FIT_TYPE'] = ptcFitType
outputStatistics[ampName]['PTC_ROW_MEAN_VARIANCE'] = inputCalib.rowMeanVariance[ampName].tolist()
outputStatistics[ampName]['PTC_MAX_RAW_MEANS'] = float(np.max(inputCalib.rawMeans[ampName]))
# To plot Covs[ij] vs flux
rawFlux = inputCalib.rawMeans[ampName].tolist()
outputStatistics[ampName]['PTC_RAW_MEANS'] = rawFlux
mask = inputCalib.expIdMask[ampName].tolist()
outputStatistics[ampName]['PTC_EXP_ID_MASK'] = mask
covs = inputCalib.covariances[ampName]
outputStatistics[ampName]['PTC_COV_10'] = covs[:, 1, 0].tolist()
outputStatistics[ampName]['PTC_COV_01'] = covs[:, 0, 1].tolist()
outputStatistics[ampName]['PTC_COV_11'] = covs[:, 1, 1].tolist()
outputStatistics[ampName]['PTC_COV_20'] = covs[:, 2, 0].tolist()
outputStatistics[ampName]['PTC_COV_02'] = covs[:, 0, 2].tolist()
# Calculate and save the slopes and offsets from Covs[ij] vs flux
keys = ['PTC_COV_10', 'PTC_COV_01', 'PTC_COV_11', 'PTC_COV_20',
'PTC_COV_02']
maskedFlux = np.array(rawFlux)[mask]
for key in keys:
maskedCov = np.array(outputStatistics[ampName][key])[mask]
linearFit = least_squares(modelResidual, [1., 0.0],
args=(np.array(maskedFlux), np.array(maskedCov)),
loss='cauchy')
slopeKey = key + '_FIT_SLOPE'
offsetKey = key + '_FIT_OFFSET'
successKey = key + '_FIT_SUCCESS'
outputStatistics[ampName][slopeKey] = float(linearFit.x[0])
outputStatistics[ampName][offsetKey] = float(linearFit.x[1])
outputStatistics[ampName][successKey] = linearFit.success

if ptcFitType == 'EXPAPPROXIMATION':
outputStatistics[ampName]['PTC_BFE_A00'] = float(inputCalib.ptcFitPars[ampName][0])
if ptcFitType == 'FULLCOVARIANCE':
Expand All @@ -149,7 +189,7 @@ def amplifierStatistics(self, inputCalib, camera=None):
slope = sum(rowMeanVar) / sum(2.*signal/numCols)
except ZeroDivisionError:
slope = np.nan
outputStatistics[ampName]['ROW_MEAN_VARIANCE_SLOPE'] = float(slope)
outputStatistics[ampName]['PTC_ROW_MEAN_VARIANCE_SLOPE'] = float(slope)

return outputStatistics

Expand Down Expand Up @@ -235,24 +275,4 @@ def verify(self, calib, statisticsDict, camera=None):

verifyStats[ampName] = verify

# Loop over amps to make a detector summary.
verifyDetStats = {'PTC_GAIN': [], 'PTC_NOISE': [], 'PTC_TURNOFF': [], 'PTC_BFE_A00': []}
for amp in verifyStats:
for testName in verifyStats[amp]:
if testName == 'SUCCESS':
continue
verifyDetStats[testName].append(verifyStats[amp][testName])

# If ptc model did not fit for a00 (e.g., POLYNOMIAL)
if not len(verifyDetStats['PTC_BFE_A00']):
verifyDetStats.pop('PTC_BFE_A00')

# VerifyDetStatsFinal has final boolean test over all amps
verifyDetStatsFinal = {}
for testName in verifyDetStats:
testBool = bool(np.all(list(verifyDetStats[testName])))
# Save the tests that failed
if not testBool:
verifyDetStatsFinal[testName] = bool(np.all(list(verifyDetStats[testName])))

return verifyDetStatsFinal, bool(success)
return {'AMP': verifyStats}, bool(success)
Loading