Coverage for python/lsst/ip/isr/linearize.py: 19%
228 statements
« prev ^ index » next coverage.py v6.4.1, created at 2022-07-09 06:37 -0700
« prev ^ index » next coverage.py v6.4.1, created at 2022-07-09 06:37 -0700
1#
2# LSST Data Management System
3# Copyright 2016 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 <http://www.lsstcorp.org/LegalNotices/>.
21#
22import abc
23import numpy as np
25from astropy.table import Table
27import lsst.afw.math as afwMath
28from lsst.pipe.base import Struct
29from lsst.geom import Box2I, Point2I, Extent2I
30from .applyLookupTable import applyLookupTable
31from .calibType import IsrCalib
33__all__ = ["Linearizer",
34 "LinearizeBase", "LinearizeLookupTable", "LinearizeSquared",
35 "LinearizeProportional", "LinearizePolynomial", "LinearizeSpline", "LinearizeNone"]
38class Linearizer(IsrCalib):
39 """Parameter set for linearization.
41 These parameters are included in cameraGeom.Amplifier, but
42 should be accessible externally to allow for testing.
44 Parameters
45 ----------
46 table : `numpy.array`, optional
47 Lookup table; a 2-dimensional array of floats:
48 - one row for each row index (value of coef[0] in the amplifier)
49 - one column for each image value
50 To avoid copying the table the last index should vary fastest
51 (numpy default "C" order)
52 detector : `lsst.afw.cameraGeom.Detector`, optional
53 Detector object. Passed to self.fromDetector() on init.
54 log : `logging.Logger`, optional
55 Logger to handle messages.
56 kwargs : `dict`, optional
57 Other keyword arguments to pass to the parent init.
59 Raises
60 ------
61 RuntimeError :
62 Raised if the supplied table is not 2D, or if the table has fewer
63 columns than rows (indicating that the indices are swapped).
65 Notes
66 -----
67 The linearizer attributes stored are:
69 hasLinearity : `bool`
70 Whether a linearity correction is defined for this detector.
71 override : `bool`
72 Whether the detector parameters should be overridden.
73 ampNames : `list` [`str`]
74 List of amplifier names to correct.
75 linearityCoeffs : `dict` [`str`, `numpy.array`]
76 Coefficients to use in correction. Indexed by amplifier
77 names. The format of the array depends on the type of
78 correction to apply.
79 linearityType : `dict` [`str`, `str`]
80 Type of correction to use, indexed by amplifier names.
81 linearityBBox : `dict` [`str`, `lsst.geom.Box2I`]
82 Bounding box the correction is valid over, indexed by
83 amplifier names.
84 fitParams : `dict` [`str`, `numpy.array`], optional
85 Linearity fit parameters used to construct the correction
86 coefficients, indexed as above.
87 fitParamsErr : `dict` [`str`, `numpy.array`], optional
88 Uncertainty values of the linearity fit parameters used to
89 construct the correction coefficients, indexed as above.
90 fitChiSq : `dict` [`str`, `float`], optional
91 Chi-squared value of the linearity fit, indexed as above.
92 fitResiduals : `dict` [`str`, `numpy.array`], optional
93 Residuals of the fit, indexed as above. Used for
94 calculating photdiode corrections
95 linearFit : The linear fit to the low flux region of the curve.
96 [intercept, slope].
97 tableData : `numpy.array`, optional
98 Lookup table data for the linearity correction.
99 """
100 _OBSTYPE = "LINEARIZER"
101 _SCHEMA = 'Gen3 Linearizer'
102 _VERSION = 1.1
104 def __init__(self, table=None, **kwargs):
105 self.hasLinearity = False
106 self.override = False
108 self.ampNames = list()
109 self.linearityCoeffs = dict()
110 self.linearityType = dict()
111 self.linearityBBox = dict()
112 self.fitParams = dict()
113 self.fitParamsErr = dict()
114 self.fitChiSq = dict()
115 self.fitResiduals = dict()
116 self.linearFit = dict()
117 self.tableData = None
118 if table is not None:
119 if len(table.shape) != 2:
120 raise RuntimeError("table shape = %s; must have two dimensions" % (table.shape,))
121 if table.shape[1] < table.shape[0]:
122 raise RuntimeError("table shape = %s; indices are switched" % (table.shape,))
123 self.tableData = np.array(table, order="C")
125 super().__init__(**kwargs)
126 self.requiredAttributes.update(['hasLinearity', 'override',
127 'ampNames',
128 'linearityCoeffs', 'linearityType', 'linearityBBox',
129 'fitParams', 'fitParamsErr', 'fitChiSq',
130 'fitResiduals', 'linearFit', 'tableData'])
132 def updateMetadata(self, setDate=False, **kwargs):
133 """Update metadata keywords with new values.
135 This calls the base class's method after ensuring the required
136 calibration keywords will be saved.
138 Parameters
139 ----------
140 setDate : `bool`, optional
141 Update the CALIBDATE fields in the metadata to the current
142 time. Defaults to False.
143 kwargs :
144 Other keyword parameters to set in the metadata.
145 """
146 kwargs['HAS_LINEARITY'] = self.hasLinearity
147 kwargs['OVERRIDE'] = self.override
148 kwargs['HAS_TABLE'] = self.tableData is not None
150 super().updateMetadata(setDate=setDate, **kwargs)
152 def fromDetector(self, detector):
153 """Read linearity parameters from a detector.
155 Parameters
156 ----------
157 detector : `lsst.afw.cameraGeom.detector`
158 Input detector with parameters to use.
160 Returns
161 -------
162 calib : `lsst.ip.isr.Linearizer`
163 The calibration constructed from the detector.
164 """
165 self._detectorName = detector.getName()
166 self._detectorSerial = detector.getSerial()
167 self._detectorId = detector.getId()
168 self.hasLinearity = True
170 # Do not translate Threshold, Maximum, Units.
171 for amp in detector.getAmplifiers():
172 ampName = amp.getName()
173 self.ampNames.append(ampName)
174 self.linearityType[ampName] = amp.getLinearityType()
175 self.linearityCoeffs[ampName] = amp.getLinearityCoeffs()
176 self.linearityBBox[ampName] = amp.getBBox()
178 return self
180 @classmethod
181 def fromDict(cls, dictionary):
182 """Construct a calibration from a dictionary of properties
184 Parameters
185 ----------
186 dictionary : `dict`
187 Dictionary of properties
189 Returns
190 -------
191 calib : `lsst.ip.isr.Linearity`
192 Constructed calibration.
194 Raises
195 ------
196 RuntimeError
197 Raised if the supplied dictionary is for a different
198 calibration.
199 """
201 calib = cls()
203 if calib._OBSTYPE != dictionary['metadata']['OBSTYPE']:
204 raise RuntimeError(f"Incorrect linearity supplied. Expected {calib._OBSTYPE}, "
205 f"found {dictionary['metadata']['OBSTYPE']}")
207 calib.setMetadata(dictionary['metadata'])
209 calib.hasLinearity = dictionary.get('hasLinearity',
210 dictionary['metadata'].get('HAS_LINEARITY', False))
211 calib.override = dictionary.get('override', True)
213 if calib.hasLinearity:
214 for ampName in dictionary['amplifiers']:
215 amp = dictionary['amplifiers'][ampName]
216 calib.ampNames.append(ampName)
217 calib.linearityCoeffs[ampName] = np.array(amp.get('linearityCoeffs', [0.0]))
218 calib.linearityType[ampName] = amp.get('linearityType', 'None')
219 calib.linearityBBox[ampName] = amp.get('linearityBBox', None)
221 calib.fitParams[ampName] = np.array(amp.get('fitParams', [0.0]))
222 calib.fitParamsErr[ampName] = np.array(amp.get('fitParamsErr', [0.0]))
223 calib.fitChiSq[ampName] = amp.get('fitChiSq', np.nan)
224 calib.fitResiduals[ampName] = np.array(amp.get('fitResiduals', [0.0]))
225 calib.linearFit[ampName] = np.array(amp.get('linearFit', [0.0]))
227 calib.tableData = dictionary.get('tableData', None)
228 if calib.tableData:
229 calib.tableData = np.array(calib.tableData)
231 return calib
233 def toDict(self):
234 """Return linearity parameters as a dict.
236 Returns
237 -------
238 outDict : `dict`:
239 """
240 self.updateMetadata()
242 outDict = {'metadata': self.getMetadata(),
243 'detectorName': self._detectorName,
244 'detectorSerial': self._detectorSerial,
245 'detectorId': self._detectorId,
246 'hasTable': self.tableData is not None,
247 'amplifiers': dict(),
248 }
249 for ampName in self.linearityType:
250 outDict['amplifiers'][ampName] = {'linearityType': self.linearityType[ampName],
251 'linearityCoeffs': self.linearityCoeffs[ampName].tolist(),
252 'linearityBBox': self.linearityBBox[ampName],
253 'fitParams': self.fitParams[ampName].tolist(),
254 'fitParamsErr': self.fitParamsErr[ampName].tolist(),
255 'fitChiSq': self.fitChiSq[ampName],
256 'fitResiduals': self.fitResiduals[ampName].tolist(),
257 'linearFit': self.linearFit[ampName].tolist()}
258 if self.tableData is not None:
259 outDict['tableData'] = self.tableData.tolist()
261 return outDict
263 @classmethod
264 def fromTable(cls, tableList):
265 """Read linearity from a FITS file.
267 This method uses the `fromDict` method to create the
268 calibration, after constructing an appropriate dictionary from
269 the input tables.
271 Parameters
272 ----------
273 tableList : `list` [`astropy.table.Table`]
274 afwTable read from input file name.
276 Returns
277 -------
278 linearity : `~lsst.ip.isr.linearize.Linearizer``
279 Linearity parameters.
281 Notes
282 -----
283 The method reads a FITS file with 1 or 2 extensions. The metadata is
284 read from the header of extension 1, which must exist. Then the table
285 is loaded, and the ['AMPLIFIER_NAME', 'TYPE', 'COEFFS', 'BBOX_X0',
286 'BBOX_Y0', 'BBOX_DX', 'BBOX_DY'] columns are read and used to set each
287 dictionary by looping over rows.
288 Extension 2 is then attempted to read in the try block (which only
289 exists for lookup tables). It has a column named 'LOOKUP_VALUES' that
290 contains a vector of the lookup entries in each row.
291 """
292 coeffTable = tableList[0]
294 metadata = coeffTable.meta
295 inDict = dict()
296 inDict['metadata'] = metadata
297 inDict['hasLinearity'] = metadata.get('HAS_LINEARITY', False)
298 inDict['amplifiers'] = dict()
300 for record in coeffTable:
301 ampName = record['AMPLIFIER_NAME']
303 fitParams = record['FIT_PARAMS'] if 'FIT_PARAMS' in record.columns else np.array([0.0])
304 fitParamsErr = record['FIT_PARAMS_ERR'] if 'FIT_PARAMS_ERR' in record.columns else np.array([0.0])
305 fitChiSq = record['RED_CHI_SQ'] if 'RED_CHI_SQ' in record.columns else np.nan
306 fitResiduals = record['FIT_RES'] if 'FIT_RES' in record.columns else np.array([0.0])
307 linearFit = record['LIN_FIT'] if 'LIN_FIT' in record.columns else np.array([0.0])
309 inDict['amplifiers'][ampName] = {
310 'linearityType': record['TYPE'],
311 'linearityCoeffs': record['COEFFS'],
312 'linearityBBox': Box2I(Point2I(record['BBOX_X0'], record['BBOX_Y0']),
313 Extent2I(record['BBOX_DX'], record['BBOX_DY'])),
314 'fitParams': fitParams,
315 'fitParamsErr': fitParamsErr,
316 'fitChiSq': fitChiSq,
317 'fitResiduals': fitResiduals,
318 'linearFit': linearFit,
319 }
321 if len(tableList) > 1:
322 tableData = tableList[1]
323 inDict['tableData'] = [record['LOOKUP_VALUES'] for record in tableData]
325 return cls().fromDict(inDict)
327 def toTable(self):
328 """Construct a list of tables containing the information in this
329 calibration.
331 The list of tables should create an identical calibration
332 after being passed to this class's fromTable method.
334 Returns
335 -------
336 tableList : `list` [`astropy.table.Table`]
337 List of tables containing the linearity calibration
338 information.
339 """
341 tableList = []
342 self.updateMetadata()
343 catalog = Table([{'AMPLIFIER_NAME': ampName,
344 'TYPE': self.linearityType[ampName],
345 'COEFFS': self.linearityCoeffs[ampName],
346 'BBOX_X0': self.linearityBBox[ampName].getMinX(),
347 'BBOX_Y0': self.linearityBBox[ampName].getMinY(),
348 'BBOX_DX': self.linearityBBox[ampName].getWidth(),
349 'BBOX_DY': self.linearityBBox[ampName].getHeight(),
350 'FIT_PARAMS': self.fitParams[ampName],
351 'FIT_PARAMS_ERR': self.fitParamsErr[ampName],
352 'RED_CHI_SQ': self.fitChiSq[ampName],
353 'FIT_RES': self.fitResiduals[ampName],
354 'LIN_FIT': self.linearFit[ampName],
355 } for ampName in self.ampNames])
356 catalog.meta = self.getMetadata().toDict()
357 tableList.append(catalog)
359 if self.tableData is not None:
360 catalog = Table([{'LOOKUP_VALUES': value} for value in self.tableData])
361 tableList.append(catalog)
362 return(tableList)
364 def getLinearityTypeByName(self, linearityTypeName):
365 """Determine the linearity class to use from the type name.
367 Parameters
368 ----------
369 linearityTypeName : str
370 String name of the linearity type that is needed.
372 Returns
373 -------
374 linearityType : `~lsst.ip.isr.linearize.LinearizeBase`
375 The appropriate linearity class to use. If no matching class
376 is found, `None` is returned.
377 """
378 for t in [LinearizeLookupTable,
379 LinearizeSquared,
380 LinearizePolynomial,
381 LinearizeProportional,
382 LinearizeSpline,
383 LinearizeNone]:
384 if t.LinearityType == linearityTypeName:
385 return t
386 return None
388 def validate(self, detector=None, amplifier=None):
389 """Validate linearity for a detector/amplifier.
391 Parameters
392 ----------
393 detector : `lsst.afw.cameraGeom.Detector`, optional
394 Detector to validate, along with its amplifiers.
395 amplifier : `lsst.afw.cameraGeom.Amplifier`, optional
396 Single amplifier to validate.
398 Raises
399 ------
400 RuntimeError :
401 Raised if there is a mismatch in linearity parameters, and
402 the cameraGeom parameters are not being overridden.
403 """
404 amplifiersToCheck = []
405 if detector:
406 if self._detectorName != detector.getName():
407 raise RuntimeError("Detector names don't match: %s != %s" %
408 (self._detectorName, detector.getName()))
409 if int(self._detectorId) != int(detector.getId()):
410 raise RuntimeError("Detector IDs don't match: %s != %s" %
411 (int(self._detectorId), int(detector.getId())))
412 if self._detectorSerial != detector.getSerial():
413 raise RuntimeError("Detector serial numbers don't match: %s != %s" %
414 (self._detectorSerial, detector.getSerial()))
415 if len(detector.getAmplifiers()) != len(self.linearityCoeffs.keys()):
416 raise RuntimeError("Detector number of amps = %s does not match saved value %s" %
417 (len(detector.getAmplifiers()),
418 len(self.linearityCoeffs.keys())))
419 amplifiersToCheck.extend(detector.getAmplifiers())
421 if amplifier:
422 amplifiersToCheck.extend(amplifier)
424 for amp in amplifiersToCheck:
425 ampName = amp.getName()
426 if ampName not in self.linearityCoeffs.keys():
427 raise RuntimeError("Amplifier %s is not in linearity data" %
428 (ampName, ))
429 if amp.getLinearityType() != self.linearityType[ampName]:
430 if self.override:
431 self.log.warning("Overriding amplifier defined linearityType (%s) for %s",
432 self.linearityType[ampName], ampName)
433 else:
434 raise RuntimeError("Amplifier %s type %s does not match saved value %s" %
435 (ampName, amp.getLinearityType(), self.linearityType[ampName]))
436 if (amp.getLinearityCoeffs().shape != self.linearityCoeffs[ampName].shape or not
437 np.allclose(amp.getLinearityCoeffs(), self.linearityCoeffs[ampName], equal_nan=True)):
438 if self.override:
439 self.log.warning("Overriding amplifier defined linearityCoeffs (%s) for %s",
440 self.linearityCoeffs[ampName], ampName)
441 else:
442 raise RuntimeError("Amplifier %s coeffs %s does not match saved value %s" %
443 (ampName, amp.getLinearityCoeffs(), self.linearityCoeffs[ampName]))
445 def applyLinearity(self, image, detector=None, log=None):
446 """Apply the linearity to an image.
448 If the linearity parameters are populated, use those,
449 otherwise use the values from the detector.
451 Parameters
452 ----------
453 image : `~lsst.afw.image.image`
454 Image to correct.
455 detector : `~lsst.afw.cameraGeom.detector`
456 Detector to use for linearity parameters if not already
457 populated.
458 log : `~logging.Logger`, optional
459 Log object to use for logging.
460 """
461 if log is None:
462 log = self.log
463 if detector and not self.hasLinearity:
464 self.fromDetector(detector)
466 self.validate(detector)
468 numAmps = 0
469 numLinearized = 0
470 numOutOfRange = 0
471 for ampName in self.linearityType.keys():
472 linearizer = self.getLinearityTypeByName(self.linearityType[ampName])
473 numAmps += 1
474 if linearizer is not None:
475 ampView = image.Factory(image, self.linearityBBox[ampName])
476 success, outOfRange = linearizer()(ampView, **{'coeffs': self.linearityCoeffs[ampName],
477 'table': self.tableData,
478 'log': self.log})
479 numOutOfRange += outOfRange
480 if success:
481 numLinearized += 1
482 elif log is not None:
483 log.warning("Amplifier %s did not linearize.",
484 ampName)
485 return Struct(
486 numAmps=numAmps,
487 numLinearized=numLinearized,
488 numOutOfRange=numOutOfRange
489 )
492class LinearizeBase(metaclass=abc.ABCMeta):
493 """Abstract base class functor for correcting non-linearity.
495 Subclasses must define __call__ and set class variable
496 LinearityType to a string that will be used for linearity type in
497 the cameraGeom.Amplifier.linearityType field.
499 All linearity corrections should be defined in terms of an
500 additive correction, such that:
502 corrected_value = uncorrected_value + f(uncorrected_value)
503 """
504 LinearityType = None # linearity type, a string used for AmpInfoCatalogs
506 @abc.abstractmethod
507 def __call__(self, image, **kwargs):
508 """Correct non-linearity.
510 Parameters
511 ----------
512 image : `lsst.afw.image.Image`
513 Image to be corrected
514 kwargs : `dict`
515 Dictionary of parameter keywords:
516 ``"coeffs"``
517 Coefficient vector (`list` or `numpy.array`).
518 ``"table"``
519 Lookup table data (`numpy.array`).
520 ``"log"``
521 Logger to handle messages (`logging.Logger`).
523 Returns
524 -------
525 output : `bool`
526 If true, a correction was applied successfully.
528 Raises
529 ------
530 RuntimeError:
531 Raised if the linearity type listed in the
532 detector does not match the class type.
533 """
534 pass
537class LinearizeLookupTable(LinearizeBase):
538 """Correct non-linearity with a persisted lookup table.
540 The lookup table consists of entries such that given
541 "coefficients" c0, c1:
543 for each i,j of image:
544 rowInd = int(c0)
545 colInd = int(c1 + uncorrImage[i,j])
546 corrImage[i,j] = uncorrImage[i,j] + table[rowInd, colInd]
548 - c0: row index; used to identify which row of the table to use
549 (typically one per amplifier, though one can have multiple
550 amplifiers use the same table)
551 - c1: column index offset; added to the uncorrected image value
552 before truncation; this supports tables that can handle
553 negative image values; also, if the c1 ends with .5 then
554 the nearest index is used instead of truncating to the
555 next smaller index
556 """
557 LinearityType = "LookupTable"
559 def __call__(self, image, **kwargs):
560 """Correct for non-linearity.
562 Parameters
563 ----------
564 image : `lsst.afw.image.Image`
565 Image to be corrected
566 kwargs : `dict`
567 Dictionary of parameter keywords:
568 ``"coeffs"``
569 Columnation vector (`list` or `numpy.array`).
570 ``"table"``
571 Lookup table data (`numpy.array`).
572 ``"log"``
573 Logger to handle messages (`logging.Logger`).
575 Returns
576 -------
577 output : `tuple` [`bool`, `int`]
578 If true, a correction was applied successfully. The
579 integer indicates the number of pixels that were
580 uncorrectable by being out of range.
582 Raises
583 ------
584 RuntimeError:
585 Raised if the requested row index is out of the table
586 bounds.
587 """
588 numOutOfRange = 0
590 rowInd, colIndOffset = kwargs['coeffs'][0:2]
591 table = kwargs['table']
592 log = kwargs['log']
594 numTableRows = table.shape[0]
595 rowInd = int(rowInd)
596 if rowInd < 0 or rowInd > numTableRows:
597 raise RuntimeError("LinearizeLookupTable rowInd=%s not in range[0, %s)" %
598 (rowInd, numTableRows))
599 tableRow = np.array(table[rowInd, :], dtype=image.getArray().dtype)
601 numOutOfRange += applyLookupTable(image, tableRow, colIndOffset)
603 if numOutOfRange > 0 and log is not None:
604 log.warning("%s pixels were out of range of the linearization table",
605 numOutOfRange)
606 if numOutOfRange < image.getArray().size:
607 return True, numOutOfRange
608 else:
609 return False, numOutOfRange
612class LinearizePolynomial(LinearizeBase):
613 """Correct non-linearity with a polynomial mode.
615 corrImage = uncorrImage + sum_i c_i uncorrImage^(2 + i)
617 where c_i are the linearity coefficients for each amplifier.
618 Lower order coefficients are not included as they duplicate other
619 calibration parameters:
620 ``"k0"``
621 A coefficient multiplied by uncorrImage**0 is equivalent to
622 bias level. Irrelevant for correcting non-linearity.
623 ``"k1"``
624 A coefficient multiplied by uncorrImage**1 is proportional
625 to the gain. Not necessary for correcting non-linearity.
626 """
627 LinearityType = "Polynomial"
629 def __call__(self, image, **kwargs):
630 """Correct non-linearity.
632 Parameters
633 ----------
634 image : `lsst.afw.image.Image`
635 Image to be corrected
636 kwargs : `dict`
637 Dictionary of parameter keywords:
638 ``"coeffs"``
639 Coefficient vector (`list` or `numpy.array`).
640 If the order of the polynomial is n, this list
641 should have a length of n-1 ("k0" and "k1" are
642 not needed for the correction).
643 ``"log"``
644 Logger to handle messages (`logging.Logger`).
646 Returns
647 -------
648 output : `tuple` [`bool`, `int`]
649 If true, a correction was applied successfully. The
650 integer indicates the number of pixels that were
651 uncorrectable by being out of range.
652 """
653 if not np.any(np.isfinite(kwargs['coeffs'])):
654 return False, 0
655 if not np.any(kwargs['coeffs']):
656 return False, 0
658 ampArray = image.getArray()
659 correction = np.zeros_like(ampArray)
660 for order, coeff in enumerate(kwargs['coeffs'], start=2):
661 correction += coeff * np.power(ampArray, order)
662 ampArray += correction
664 return True, 0
667class LinearizeSquared(LinearizeBase):
668 """Correct non-linearity with a squared model.
670 corrImage = uncorrImage + c0*uncorrImage^2
672 where c0 is linearity coefficient 0 for each amplifier.
673 """
674 LinearityType = "Squared"
676 def __call__(self, image, **kwargs):
677 """Correct for non-linearity.
679 Parameters
680 ----------
681 image : `lsst.afw.image.Image`
682 Image to be corrected
683 kwargs : `dict`
684 Dictionary of parameter keywords:
685 ``"coeffs"``
686 Coefficient vector (`list` or `numpy.array`).
687 ``"log"``
688 Logger to handle messages (`logging.Logger`).
690 Returns
691 -------
692 output : `tuple` [`bool`, `int`]
693 If true, a correction was applied successfully. The
694 integer indicates the number of pixels that were
695 uncorrectable by being out of range.
696 """
698 sqCoeff = kwargs['coeffs'][0]
699 if sqCoeff != 0:
700 ampArr = image.getArray()
701 ampArr *= (1 + sqCoeff*ampArr)
702 return True, 0
703 else:
704 return False, 0
707class LinearizeSpline(LinearizeBase):
708 """Correct non-linearity with a spline model.
710 corrImage = uncorrImage - Spline(coeffs, uncorrImage)
712 Notes
713 -----
715 The spline fit calculates a correction as a function of the
716 expected linear flux term. Because of this, the correction needs
717 to be subtracted from the observed flux.
719 """
720 LinearityType = "Spline"
722 def __call__(self, image, **kwargs):
723 """Correct for non-linearity.
725 Parameters
726 ----------
727 image : `lsst.afw.image.Image`
728 Image to be corrected
729 kwargs : `dict`
730 Dictionary of parameter keywords:
731 ``"coeffs"``
732 Coefficient vector (`list` or `numpy.array`).
733 ``"log"``
734 Logger to handle messages (`logging.Logger`).
736 Returns
737 -------
738 output : `tuple` [`bool`, `int`]
739 If true, a correction was applied successfully. The
740 integer indicates the number of pixels that were
741 uncorrectable by being out of range.
742 """
743 splineCoeff = kwargs['coeffs']
744 centers, values = np.split(splineCoeff, 2)
745 interp = afwMath.makeInterpolate(centers.tolist(), values.tolist(),
746 afwMath.stringToInterpStyle("AKIMA_SPLINE"))
748 ampArr = image.getArray()
749 delta = interp.interpolate(ampArr.flatten())
750 ampArr -= np.array(delta).reshape(ampArr.shape)
752 return True, 0
755class LinearizeProportional(LinearizeBase):
756 """Do not correct non-linearity.
757 """
758 LinearityType = "Proportional"
760 def __call__(self, image, **kwargs):
761 """Do not correct for non-linearity.
763 Parameters
764 ----------
765 image : `lsst.afw.image.Image`
766 Image to be corrected
767 kwargs : `dict`
768 Dictionary of parameter keywords:
769 ``"coeffs"``
770 Coefficient vector (`list` or `numpy.array`).
771 ``"log"``
772 Logger to handle messages (`logging.Logger`).
774 Returns
775 -------
776 output : `tuple` [`bool`, `int`]
777 If true, a correction was applied successfully. The
778 integer indicates the number of pixels that were
779 uncorrectable by being out of range.
780 """
781 return True, 0
784class LinearizeNone(LinearizeBase):
785 """Do not correct non-linearity.
786 """
787 LinearityType = "None"
789 def __call__(self, image, **kwargs):
790 """Do not correct for non-linearity.
792 Parameters
793 ----------
794 image : `lsst.afw.image.Image`
795 Image to be corrected
796 kwargs : `dict`
797 Dictionary of parameter keywords:
798 ``"coeffs"``
799 Coefficient vector (`list` or `numpy.array`).
800 ``"log"``
801 Logger to handle messages (`logging.Logger`).
803 Returns
804 -------
805 output : `tuple` [`bool`, `int`]
806 If true, a correction was applied successfully. The
807 integer indicates the number of pixels that were
808 uncorrectable by being out of range.
809 """
810 return True, 0