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