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-47263: Simplify amp offset correction logging and config #358

Open
wants to merge 1 commit into
base: main
Choose a base branch
from
Open
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
136 changes: 100 additions & 36 deletions python/lsst/ip/isr/ampOffset.py
Original file line number Diff line number Diff line change
Expand Up @@ -214,10 +214,10 @@ def run(self, exposure):

# Safety check: do any pixels remain for amp offset estimation?
if (exp.mask.array & bitMask).all():
log_fn = self.log.warning if self.config.doApplyAmpOffset else self.log.debug
log_fn = self.log.warning if self.config.doApplyAmpOffset else self.log.info
log_fn(
"All pixels masked: cannot calculate any amp offset corrections. All pedestals are being set "
"to zero."
"All pixels masked: cannot calculate any amp offset corrections. "
"All pedestals are being set to zero."
)
pedestals = np.zeros(len(amps))
else:
Expand All @@ -243,8 +243,8 @@ def run(self, exposure):

metadata = exposure.getMetadata() # Exposure metadata.
self.metadata["AMPOFFSET_PEDESTALS"] = {} # Task metadata.
for amp, pedestal in zip(amps, pedestals):
ampName = amp.getName()
ampNames = [amp.getName() for amp in amps]
for ampName, amp, pedestal in zip(ampNames, amps, pedestals):
# Add the amp pedestal to the exposure metadata.
metadata.set(
f"LSST ISR AMPOFFSET PEDESTAL {ampName}",
Expand All @@ -259,15 +259,16 @@ def run(self, exposure):
self.metadata["AMPOFFSET_PEDESTALS"][ampName] = float(pedestal)
if self.config.doApplyAmpOffset:
status = "subtracted from exposure"
metadata.set(
"LSST ISR AMPOFFSET PEDESTAL SUBTRACTED", True, "Amp pedestals have been subtracted"
)
metadata.set("LSST ISR AMPOFFSET PEDESTAL SUBTRACTED", True, "Amp pedestals have been subtracted")
else:
status = "not subtracted from exposure"
metadata.set(
"LSST ISR AMPOFFSET PEDESTAL SUBTRACTED", False, "Amp pedestals have not been subtracted"
)
self.log.info(f"amp pedestal values ({status}): {', '.join([f'{x:.4f}' for x in pedestals])}")
ampPedestalReport = ", ".join(
[f"{ampName}: {ampPedestal:.4f}" for (ampName, ampPedestal) in zip(ampNames, pedestals)]
)
self.log.info(f"Amp pedestal values ({status}): {ampPedestalReport}")

return Struct(pedestals=pedestals)

Expand Down Expand Up @@ -344,8 +345,9 @@ def getAmpAssociations(self, amps):
if not np.all(ampAssociations == ampAssociations.T):
raise RuntimeError("The `ampAssociations` is not symmetric about the diagonal.")

self.log.debug("amp associations:\n%s", ampAssociations)
self.log.debug("amp sides:\n%s", ampSides)
with np.printoptions(linewidth=200):
self.log.debug("amp associations:\n%s", ampAssociations)
self.log.debug("amp sides:\n%s", ampSides)

return ampAssociations, ampSides

Expand Down Expand Up @@ -408,8 +410,14 @@ def getAmpOffsets(self, im, amps, associations, sides):
"""
ampsOffsets = np.zeros(len(amps))
ampsEdges = self.getAmpEdges(im, amps, sides)
ampsNames = [amp.getName() for amp in amps]
interfaceOffsetLookup = {}

interfaceIds = []
interfaceOffsetOriginals = []
ampEdgeGoodFracs = []
minFracFails = []
maxOffsetFails = []
for ampId, ampAssociations in enumerate(associations):
ampNeighbors = np.ravel(np.where(ampAssociations < 0))
for ampNeighbor in ampNeighbors:
Expand All @@ -420,11 +428,60 @@ def getAmpOffsets(self, im, amps, associations, sides):
edgeA = ampsEdges[ampId][ampSide]
edgeB = ampsEdges[ampNeighbor][(ampSide + 2) % 4]
if ampId < ampNeighbor:
interfaceOffset = self.getInterfaceOffset(ampId, ampNeighbor, edgeA, edgeB)
interfaceOffsetLookup[f"{ampId}{ampNeighbor}"] = interfaceOffset
(
interfaceId,
interfaceOffset,
interfaceOffsetOriginal,
ampEdgeGoodFrac,
minFracFail,
maxOffsetFail,
) = self.getInterfaceOffset(ampsNames[ampId], ampsNames[ampNeighbor], edgeA, edgeB)
interfaceIds.append(interfaceId)
interfaceOffsetOriginals.append(interfaceOffsetOriginal)
ampEdgeGoodFracs.append(ampEdgeGoodFrac)
minFracFails.append(minFracFail)
maxOffsetFails.append(maxOffsetFail)
interfaceOffsetLookup[f"{ampId:02d}:{ampNeighbor:02d}"] = interfaceOffset
else:
interfaceOffset = -interfaceOffsetLookup[f"{ampNeighbor}{ampId}"]
interfaceOffset = -interfaceOffsetLookup[f"{ampNeighbor:02d}:{ampId:02d}"]
ampsOffsets[ampId] += interfaceWeight * interfaceOffset
if interfaceOffsetOriginals:
quartile_summary = np.nanpercentile(interfaceOffsetOriginals, [0, 25, 50, 75, 100])
self.log.info(
"Amp offset quartile summary (min, Q1, Q2, Q3, max): %.4f, %.4f, %.4f, %.4f, %.4f",
*quartile_summary,
)
log_fn = self.log.warning if self.config.doApplyAmpOffset else self.log.info
if any(minFracFails):
log_fn(
"The fraction of unmasked edge pixels for the following amp interfaces is below the "
"configured threshold (%s): %s",
self.config.ampEdgeMinFrac,
", ".join(
[
f"{interfaceId} ({ampEdgeGoodFrac:0.2f})"
for interfaceId, ampEdgeGoodFrac, minFracFail in zip(
interfaceIds, ampEdgeGoodFracs, minFracFails
)
if minFracFail
]
),
)
if any(maxOffsetFails):
log_fn(
"Absolute amp offsets exceed the configured maximum (%s) and have been set to zero for the "
"following amp interfaces: %s",
self.config.ampEdgeMaxOffset,
", ".join(
[
f"{interfaceId}={np.abs(interfaceOffset):0.2f}"
for interfaceId, interfaceOffset, maxOffsetFail in zip(
interfaceIds, interfaceOffsetOriginals, maxOffsetFails
)
if maxOffsetFail
]
),
)
return ampsOffsets

def getAmpEdges(self, im, amps, ampSides):
Expand Down Expand Up @@ -474,16 +531,16 @@ def getAmpEdges(self, im, amps, ampSides):
ampEdges[ampId][ampSide] = np.nanmedian(strip, axis=ampSide % 2) # 1D medianified strip
return ampEdges

def getInterfaceOffset(self, ampIdA, ampIdB, edgeA, edgeB):
def getInterfaceOffset(self, ampNameA, ampNameB, edgeA, edgeB):
"""Calculate the amp offset for a given interface between two
amplifiers.

Parameters
----------
ampIdA : int
ID of the first amplifier.
ampIdB : int
ID of the second amplifier.
ampNameA : str
Name of the first amplifier.
ampNameB : str
Name of the second amplifier.
edgeA : numpy.ndarray
Amp edge for the first amplifier.
edgeB : numpy.ndarray
Expand All @@ -494,8 +551,19 @@ def getInterfaceOffset(self, ampIdA, ampIdB, edgeA, edgeB):
interfaceOffset : float
The calculated amp offset value for the given interface between
amps A and B.
interfaceOffsetOriginal : float
The original calculated amp offset value for the given interface
between amps A and B.
ampEdgeGoodFrac : float
Fraction of viable pixel rows along the amp edge.
minFracFail : bool
True if the fraction of unmasked pixel rows is below the
ampEdgeMinFrac threshold.
maxOffsetFail : bool
True if the absolute offset value exceeds the ampEdgeMaxOffset
threshold.
"""
interfaceId = f"{ampIdA}{ampIdB}"
interfaceId = f"{ampNameA}:{ampNameB}"
sctrl = StatisticsControl()
# NOTE: Taking the difference with the order below fixes the sign flip
# in the B matrix.
Expand All @@ -512,30 +580,26 @@ def getInterfaceOffset(self, ampIdA, ampIdB, edgeA, edgeB):
edgeDiffAvg[np.isnan(edgeDiff)] = np.nan
# Take clipped mean of rolling average data as amp offset value.
interfaceOffset = makeStatistics(edgeDiffAvg, MEANCLIP, sctrl).getValue()
interfaceOffsetOriginal = interfaceOffset
ampEdgeGoodFrac = 1 - (np.sum(np.isnan(edgeDiffAvg)) / len(edgeDiffAvg))

# Perform a couple of do-no-harm safety checks:
# a) The fraction of unmasked pixel rows is > ampEdgeMinFrac,
# b) The absolute offset ADU value is < ampEdgeMaxOffset.
minFracFail = ampEdgeGoodFrac < self.config.ampEdgeMinFrac
maxOffsetFail = np.abs(interfaceOffset) > self.config.ampEdgeMaxOffset
if minFracFail or maxOffsetFail:
log_fn = self.log.warning if self.config.doApplyAmpOffset else self.log.debug
if minFracFail:
log_fn(
f"The fraction of unmasked pixels for amp interface {interfaceId} is below the "
f"threshold ({ampEdgeGoodFrac:.2f} < {self.config.ampEdgeMinFrac}). Resetting the "
f"interface offset from {interfaceOffset} to 0."
)
if maxOffsetFail:
log_fn(
"The absolute offset value exceeds the limit "
f"({np.abs(interfaceOffset):.2f} > {self.config.ampEdgeMaxOffset} ADU). Resetting "
f"the interface offset from {interfaceOffset} to 0."
)
interfaceOffset = 0
self.log.debug(
f"amp interface {interfaceId} : "
f"viable edge difference frac = {ampEdgeGoodFrac}, "
f"interface offset = {interfaceOffset:.3f}"
f"viable edge difference frac = {ampEdgeGoodFrac:.2f}, "
f"amp offset = {interfaceOffsetOriginal:.3f}"
)
return (
interfaceId,
interfaceOffset,
interfaceOffsetOriginal,
ampEdgeGoodFrac,
minFracFail,
maxOffsetFail,
)
return interfaceOffset
4 changes: 2 additions & 2 deletions python/lsst/ip/isr/isrTask.py
Original file line number Diff line number Diff line change
Expand Up @@ -1731,9 +1731,9 @@ def run(self, ccdExposure, *, camera=None, bias=None, linearizer=None,
# Calculate amp offset corrections within the CCD.
if self.config.doAmpOffset:
if self.config.ampOffset.doApplyAmpOffset:
self.log.info("Calculating and applying amp offset corrections.")
self.log.info("Measuring amp offset corrections only, without applying them.")
else:
self.log.info("Calculating amp offset corrections without applying them.")
self.log.info("Measuring and applying amp offset corrections.")
self.ampOffset.run(ccdExposure)

if self.config.doMeasureBackground:
Expand Down
4 changes: 2 additions & 2 deletions python/lsst/ip/isr/isrTaskLSST.py
Original file line number Diff line number Diff line change
Expand Up @@ -1826,9 +1826,9 @@ def run(self, ccdExposure, *, dnlLUT=None, bias=None, deferredChargeCalib=None,
# Calculate amp offset corrections within the CCD.
if self.config.doAmpOffset:
if self.config.ampOffset.doApplyAmpOffset:
self.log.info("Calculating and applying amp offset corrections.")
self.log.info("Measuring amp offset corrections only, without applying them.")
else:
self.log.info("Calculating amp offset corrections without applying them.")
self.log.info("Measuring and applying amp offset corrections.")
self.ampOffset.run(ccdExposure)

# Calculate standard image quality statistics
Expand Down
Loading