Skip to content

Commit

Permalink
Do correct measureCti if exposure is trimmed
Browse files Browse the repository at this point in the history
  • Loading branch information
Alex-Broughton committed Oct 3, 2024
1 parent ce970b3 commit 5ed87b5
Show file tree
Hide file tree
Showing 7 changed files with 75 additions and 61 deletions.
61 changes: 18 additions & 43 deletions python/lsst/ip/isr/deferredCharge.py
Original file line number Diff line number Diff line change
Expand Up @@ -471,7 +471,7 @@ def model_results(params, signal, num_transfers, amp, start=1, stop=10, trap_typ
stop += 1

# Electronics effect optimization
output_amplifier = FloatingOutputAmplifier(v['driftscale'], v['decaytime'])
output_amplifier = FloatingOutputAmplifier(1.0, v['driftscale'], v['decaytime'])

# CTI optimization
v['cti'] = 10**v['ctiexp']
Expand Down Expand Up @@ -655,6 +655,8 @@ class FloatingOutputAmplifier:
Parameters
----------
gain : `float`
Gain of the amplifier. Currently not used.
scale : `float`
Drift scale for the amplifier.
decay_time : `float`
Expand All @@ -665,8 +667,9 @@ class FloatingOutputAmplifier:
Global CTI offset.
"""

def __init__(self, scale, decay_time, noise=0.0, offset=0.0):
def __init__(self, gain, scale, decay_time, noise=0.0, offset=0.0):

self.gain = gain
self.noise = noise
self.global_offset = offset

Expand Down Expand Up @@ -948,33 +951,13 @@ def fromTable(cls, tableList):

inDict['serialTraps'][amp] = ampTrap

# Version check
if calibVersion < 1.1:
oldCalibUnits = inDict['metadata']["UNITS"]
cls().log.warning(f"Using old version of CTI calibration ({calibVersion}), which "
f"reports units of {oldCalibUnits}, however these units are "
"likely misleading may be correct or incorrect. This version "
"also does not contain enough information to reconcile it "
"into the correct units (electron). Must assume it is in "
"units of electron and move forward.")
# Need to standardize the header units.
if oldCalibUnits == "electrons":
# The calibration *might* be in the correct units,
# but we cannot be sure.
inDict['metadata']['UNITS'] = "electron"
elif oldCalibUnits == "ADU":
# The calibration was *likely* made from images in units of adu
# Attempt to convert, but note that there is not enough
# information in the header to ensure this. The old calib
# dictionary also did not store gain information to even
# attempt to convert the only parameter (trap size) that is not
# unitless. If this calibration is used to correct an image
# during ISR, the ISR task will issue a warning and attempt to
# correct this using the gains supplied during ISR.
inDict['metadata']['UNITS'] = "adu"
else:
cls().log.warning("Reconciling units of old CTI calibration (Version "
f"{calibVersion}), but do not recognize units of "
f"\"{oldCalibUnits}\".")
#This version might be in the wrong units (not electron), and does not
# contain the gain information to convert into a a new calibration
# version.
raise RuntimeError(f"Using old version of CTI calibration (ver. {calibVersion} < 1.1), "
"which is no longer supported.")

return cls.fromDict(inDict)

Expand Down Expand Up @@ -1067,6 +1050,12 @@ class DeferredChargeConfig(Config):
doc="Number of prior pixels to use for trap correction.",
default=6,
)
useGains = Field(
dtype=bool,
doc="If true, scale by the gain.",
default=False,
deprecated="This field is no longer used. Will be removed after v28.",
)
zeroUnusedPixels = Field(
dtype=bool,
doc="If true, set serial prescan and parallel overscan to zero before correction.",
Expand Down Expand Up @@ -1119,20 +1108,6 @@ def run(self, exposure, ctiCalib, gains=None):

# Get the image and overscan units.
imageUnits = exposure.getMetadata().get("LSST ISR UNITS")
calibUnits = ctiCalib.getMetadata().get("UNITS")
calibVersion = ctiCalib.getMetadata().get("CTI_VERSION")

if calibUnits == "adu":
self.log.warning("CTI calibration is in adu, but the deferred charge correction "
"expects electron units. The calibration is likely old (version "
f"{calibVersion} < 1.1). Attempting to reconcile units, but "
"cannot guarantee that everything is correct.")
# Luckily, the trap size is the only parameter that is
# not unitless.
ctiCalib = ctiCalib.copy()
for ampName in ctiCalib.serialTraps:
ctiCalib.serialTraps[ampName].size *= gains[ampName]
ctiCalib.setMetadata(UNITS="electron")

# The deferred charge correction assumes that everything is in
# electron units. Make it so:
Expand All @@ -1153,7 +1128,7 @@ def run(self, exposure, ctiCalib, gains=None):
for amp in detector.getAmplifiers():
ampName = amp.getName()

ampImage = image[amp.getRawBBox()]
ampImage = image[amp.getRawDataBBox()]
if self.config.zeroUnusedPixels:
# We don't apply overscan subtraction, so zero these
# out for now.
Expand Down
36 changes: 36 additions & 0 deletions python/lsst/ip/isr/isrFunctions.py
Original file line number Diff line number Diff line change
Expand Up @@ -1288,3 +1288,39 @@ def getExposureReadNoises(exposure):
else:
readnoises[ampName] = amp.getReadNoise()
return readnoises


def isTrimmedExposure(exposure):
"""Check if the unused pixels (pre-/over-scan pixels) have
been trimmed from an exposure.
Parameters
----------
exposure : `lsst.afw.image.Exposure`
The exposure to check.
Returns
-------
result : `bool`
True if the image is trimmed, else False.
"""
return exposure.getDetector().getBBox() == exposure.getBBox()


def isTrimmedImage(image, detector):
"""Check if the unused pixels (pre-/over-scan pixels) have
been trimmed from an image
Parameters
----------
image : `lsst.afw.image.Image`
The image to check.
detector : `lsst.afw.cameraGeom.Detector`
The detector associated with the image.
Returns
-------
result : `bool`
True if the image is trimmed, else False.
"""
return detector.getBBox() == image.getBBox()
1 change: 1 addition & 0 deletions python/lsst/ip/isr/isrMockLSST.py
Original file line number Diff line number Diff line change
Expand Up @@ -926,6 +926,7 @@ def amplifierAddDeferredCharge(self, exposure, amp):
# Create a fake amplifier object that contains some deferred charge
# paramters.
floatingOutputAmplifier = FloatingOutputAmplifier(
gain=1.0, # Everyhting is already in electrons.
scale=driftScale,
decay_time=decayTime,
noise=0.0,
Expand Down
16 changes: 11 additions & 5 deletions python/lsst/ip/isr/isrStatistics.py
Original file line number Diff line number Diff line change
Expand Up @@ -33,7 +33,7 @@
import lsst.pex.config as pexConfig

from lsst.afw.cameraGeom import ReadoutCorner
from .isrFunctions import gainContext
from .isrFunctions import gainContext, isTrimmedExposure


class IsrStatisticsTaskConfig(pexConfig.Config):
Expand Down Expand Up @@ -324,7 +324,8 @@ def measureCti(self, inputExp, overscans, gains):
inputExp : `lsst.afw.image.Exposure`
Exposure to measure.
overscans : `list` [`lsst.pipe.base.Struct`]
List of overscan results. Expected fields are:
List of overscan results (expects base units of adu).
Expected fields are:
``imageFit``
Value or fit subtracted from the amplifier image data
Expand Down Expand Up @@ -358,16 +359,19 @@ def measureCti(self, inputExp, overscans, gains):
if imageUnits == "adu":
applyGain = True

# Check if the image is trimmed.
isTrimmed = isTrimmedExposure(inputExp)

# Ensure we have the same number of overscans as amplifiers.
assert len(overscans) == len(detector.getAmplifiers())

with gainContext(inputExp, image, applyGain, gains, isTrimmed=False):
with gainContext(inputExp, image, applyGain, gains, isTrimmed=isTrimmed:
for ampIter, amp in enumerate(detector.getAmplifiers()):
ampStats = {}
readoutCorner = amp.getReadoutCorner()

# Full data region.
dataRegion = image[amp.getRawDataBBox()]
dataRegion = image[amp.getBBox() if isTrimmed else amp.getRawDataBBox()]

# Get the mean of the image
ampStats["IMAGE_MEAN"] = afwMath.makeStatistics(dataRegion, self.statType,
Expand Down Expand Up @@ -415,7 +419,9 @@ def measureCti(self, inputExp, overscans, gains):
osMean = afwMath.makeStatistics(overscanImage.image.array[:nRows, column],
self.statType, self.statControl).getValue()
columns.append(column)
values.append(osMean)
# The overscan input is always in adu, but it only
# makes sense to measure CTI in electron units.
values.append(osMean * gains[amp.getName()])

# We want these relative to the readout corner. If that's
# on the right side, we need to swap them.
Expand Down
2 changes: 1 addition & 1 deletion python/lsst/ip/isr/isrTask.py
Original file line number Diff line number Diff line change
Expand Up @@ -1448,7 +1448,6 @@ def run(self, ccdExposure, *, camera=None, bias=None, linearizer=None,
residMean = overscanResults.residualMean[0]
residSigma = overscanResults.residualSigma[0]

self.metadata[f"LSST ISR OVERSCAN UNITS {amp.getName()}"] = "adu"
self.metadata[f"FIT MEDIAN {amp.getName()}"] = mean
self.metadata[f"FIT STDEV {amp.getName()}"] = sigma
self.log.debug(" Overscan stats for amplifer %s: %f +/- %f",
Expand Down Expand Up @@ -2269,6 +2268,7 @@ def overscanCorrection(self, ccdExposure, amp):
ampName = amp.getName()

keyBase = "LSST ISR OVERSCAN"
metadata[f"{keyBase} UNITS {amp.getName()}"] = "adu"
# Updated quantities
if isinstance(overscanResults.overscanMean, float):
# Serial overscan correction only:
Expand Down
11 changes: 5 additions & 6 deletions python/lsst/ip/isr/isrTaskLSST.py
Original file line number Diff line number Diff line change
Expand Up @@ -27,6 +27,7 @@
from .isrStatistics import IsrStatisticsTask
from .isr import maskNans
from .ptcDataset import PhotonTransferCurveDataset
from .isrFunctions import isTrimmedExposure


class IsrTaskLSSTConnections(pipeBase.PipelineTaskConnections,
Expand Down Expand Up @@ -833,10 +834,8 @@ def overscanCorrection(self, mode, detectorConfig, detector, badAmpDict, ccdExpo
results = None
else:
# This check is to confirm that we are not trying to run
# overscan on an already trimmed image. Therefore, always
# checking just the horizontal overscan bounding box is
# sufficient.
if amp.getRawHorizontalOverscanBBox().isEmpty():
# overscan on an already trimmed image.
if isTrimmedExposure(ccdExposure):
self.log.warning(
"ISR_OSCAN: No overscan region for amp %s. Not performing overscan correction.",
ampName,
Expand Down Expand Up @@ -950,8 +949,8 @@ def addVariancePlane(self, exposure, detector):
Detector with geometry info.
"""
# NOTE: this will fail if the exposure is not trimmed.
# I am not sure if this is something we need to check for
# (or how to do it efficiently).
if not isTrimmedExposure(exposure):
raise RuntimeError("Exposure must be trimmed to add variance plane.")

isElectrons = (exposure.metadata["LSST ISR UNITS"] == "electron")

Expand Down
9 changes: 3 additions & 6 deletions python/lsst/ip/isr/linearize.py
Original file line number Diff line number Diff line change
Expand Up @@ -34,6 +34,7 @@
from lsst.geom import Box2I, Point2I, Extent2I
from .applyLookupTable import applyLookupTable
from .calibType import IsrCalib
from .isrFunctions import isTrimmedImage


class Linearizer(IsrCalib):
Expand Down Expand Up @@ -496,12 +497,8 @@ def applyLinearity(self, image, detector=None, log=None, gains=None):

self.validate(detector)

isTrimmed = None
if detector:
if detector.getBBox() == image.getBBox():
isTrimmed = True
else:
isTrimmed = False
# Check if the image is trimmed.
isTrimmed = isTrimmedImage(image, detector)

numAmps = 0
numLinearized = 0
Expand Down

0 comments on commit 5ed87b5

Please sign in to comment.