Coverage for python / lsst / cp / verify / verifyBfk.py: 20%
106 statements
« prev ^ index » next coverage.py v7.13.5, created at 2026-04-22 09:10 +0000
1# This file is part of cp_verify.
2#
3# Developed for the LSST Data Management System.
4# This product includes software developed by the LSST Project
5# (http://www.lsst.org).
6# See the COPYRIGHT file at the top-level directory of this distribution
7# for details of code ownership.
8#
9# This program is free software: you can redistribute it and/or modify
10# it under the terms of the GNU General Public License as published by
11# the Free Software Foundation, either version 3 of the License, or
12# (at your option) any later version.
13#
14# This program is distributed in the hope that it will be useful,
15# but WITHOUT ANY WARRANTY; without even the implied warranty of
16# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17# GNU General Public License for more details.
18#
19# You should have received a copy of the GNU General Public License
20# along with this program. If not, see <http://www.gnu.org/licenses/>.
21import numpy as np
22import lsst.afw.table as afwTable
23import lsst.pex.config as pexConfig
24import lsst.pipe.base.connectionTypes as cT
26from scipy.optimize import least_squares
27from lsstDebug import getDebugFrame
28from .verifyStats import CpVerifyStatsConfig, CpVerifyStatsTask, CpVerifyStatsConnections
31__all__ = ['CpVerifyBfkConfig', 'CpVerifyBfkTask']
34class CpVerifyBfkConnections(CpVerifyStatsConnections,
35 dimensions={'instrument', 'visit', 'detector'}):
36 inputExp = cT.Input(
37 name="icExp",
38 doc="Input exposure to calculate statistics for.",
39 storageClass="ExposureF",
40 dimensions=["instrument", "visit", "detector"],
41 )
42 inputCatalog = cT.Input(
43 name="icSrc",
44 doc="Input catalog to calculate statistics from.",
45 storageClass="SourceCatalog",
46 dimensions=["instrument", "visit", "detector"],
47 )
48 uncorrectedCatalog = cT.Input(
49 name="uncorrectedSrc",
50 doc="Input catalog without correction applied.",
51 storageClass="SourceCatalog",
52 dimensions=["instrument", "visit", "detector"],
53 )
54 camera = cT.PrerequisiteInput(
55 name="camera",
56 storageClass="Camera",
57 doc="Input camera.",
58 dimensions=["instrument", ],
59 isCalibration=True,
60 )
62 outputStats = cT.Output(
63 name="detectorStats",
64 doc="Output statistics from cp_verify.",
65 storageClass="StructuredDataDict",
66 dimensions=["instrument", "visit", "detector"],
67 )
70class CpVerifyBfkConfig(CpVerifyStatsConfig,
71 pipelineConnections=CpVerifyBfkConnections):
72 """Inherits from base CpVerifyStatsConfig.
73 """
75 matchRadiusPix = pexConfig.Field(
76 dtype=float,
77 default=3,
78 doc=("Match radius for matching icSourceCat objects to sourceCat objects (pixels)"),
79 )
81 def setDefaults(self):
82 super().setDefaults()
83 self.stageName = 'BFK'
84 self.catalogStatKeywords = {'BRIGHT_SLOPE': None,
85 'NUM_MATCHES': None,
86 } # noqa F841
89def exponentialModel(x, A, alpha, x0):
90 """An exponential model.
91 """
92 return A * np.exp(alpha * (x - x0))
95def modelResidual(p, x, y):
96 """Model residual for fit below.
97 """
98 return y - exponentialModel(x, *p)
101class CpVerifyBfkTask(CpVerifyStatsTask):
102 """Bfk verification sub-class, implementing the verify method.
103 """
105 ConfigClass = CpVerifyBfkConfig
106 _DefaultName = 'cpVerifyBfk'
108 def catalogStatistics(self, exposure, catalog, uncorrectedCatalog, statControl):
109 """Measure the catalog statistics.
111 Parameters
112 ----------
113 exposure : `lsst.afw.image.Exposure`
114 The exposure to measure.
115 catalog : `lsst.afw.table.Table`
116 The catalog to measure.
117 uncorrectedCatalog : `lsst.afw.table.Table`
118 The uncorrected catalog to measure.
119 statControl : `lsst.afw.math.StatisticsControl`
120 Statistics control object with parameters defined by
121 the config.
123 Returns
124 -------
125 outputStatistics : `dict` [`str`, `dict` [`str`, scalar]]
126 A dictionary indexed by the amplifier name, containing
127 dictionaries of the statistics measured and their values.
128 """
129 outputStatistics = {}
131 matches = afwTable.matchXy(catalog, uncorrectedCatalog, self.config.matchRadiusPix)
132 outputStatistics['NUM_MATCHES'] = len(matches)
134 magnitude = []
135 sizeDiff = []
136 if not matches:
137 outputStatistics['MAGNITUDES'] = magnitude
138 outputStatistics['SIZE_DIFF'] = sizeDiff
139 return outputStatistics
141 xxKey = 'ext_shapeHSM_HsmSourceMoments_xx'
142 yyKey = 'ext_shapeHSM_HsmSourceMoments_yy'
143 for source, uncorrectedSource, d in matches:
144 # This uses the simple difference in source moments.
145 sourceMagnitude = -2.5 * np.log10(source.getPsfInstFlux())
146 sourceSize = source[xxKey] + source[yyKey]
147 uncorrectedSize = uncorrectedSource[xxKey] + uncorrectedSource[yyKey]
149 magnitude.append(float(sourceMagnitude))
150 sizeDiff.append(float(uncorrectedSize - sourceSize))
152 mask = np.isfinite(magnitude) * np.isfinite(sizeDiff)
153 if 'BRIGHT_SLOPE' in self.config.catalogStatKeywords:
154 exponentialFit = least_squares(modelResidual, [0.0, -0.01, 0.0],
155 args=(np.array(magnitude)[mask], np.array(sizeDiff)[mask]),
156 loss='cauchy')
158 outputStatistics['BRIGHT_SLOPE'] = float(exponentialFit.x[1])
159 outputStatistics['FIT_SUCCESS'] = exponentialFit.success
160 self.debugFit('brightSlope', magnitude, sizeDiff, exponentialFit.x)
162 outputStatistics['MAGNITUDES'] = magnitude
163 outputStatistics['SIZE_DIFF'] = sizeDiff
165 return outputStatistics
167 def debugFit(self, stepname, xVector, yVector, yModel):
168 """Debug method for linearity fitting.
169 Parameters
170 ----------
171 stepname : `str`
172 A label to use to check if we care to debug at a given
173 line of code.
174 xVector : `numpy.array`, (N,)
175 The values to use as the independent variable in the
176 linearity fit.
177 yVector : `numpy.array`, (N,)
178 The values to use as the dependent variable in the
179 linearity fit.
180 yModel : `numpy.array`, (N,)
181 The values to use as the linearized result.
182 """
183 frame = getDebugFrame(self._display, stepname)
184 if frame:
185 import matplotlib.pyplot as plt
186 figure, axes = plt.subplots(1)
188 axes.scatter(xVector, yVector, c='blue', marker='+')
189 modelX = np.arange(np.min(xVector) - 1, np.max(xVector) + 1, 0.1)
190 axes.plot(modelX, exponentialModel(modelX, *yModel), 'r-')
191 plt.xlabel("Instrumental PSF magnitude")
192 plt.ylabel("Source size trace")
193 plt.title(f"BFK slope: {yModel[0]:.3f} + {yModel[1]:.3f} m")
194 figure.show()
195 prompt = "Press Enter or c to continue [chp]..."
197 while True:
198 ans = input(prompt).lower()
199 if ans in ("", " ", "c",):
200 break
201 elif ans in ("p", ):
202 import pdb
203 pdb.set_trace()
204 elif ans in ("h", ):
205 print("[h]elp [c]ontinue [p]db e[x]itDebug")
206 plt.close()
208 def verify(self, exposure, statisticsDict):
209 """Verify that the measured statistics meet the verification criteria.
211 Parameters
212 ----------
213 exposure : `lsst.afw.image.Exposure`
214 The exposure the statistics are from.
215 statisticsDictionary : `dict` [`str`, `dict` [`str`, scalar]],
216 Dictionary of measured statistics. The inner dictionary
217 should have keys that are statistic names (`str`) with
218 values that are some sort of scalar (`int` or `float` are
219 the mostly likely types).
221 Returns
222 -------
223 outputStatistics : `dict` [`str`, `dict` [`str`, `bool`]]
224 A dictionary indexed by the amplifier name, containing
225 dictionaries of the verification criteria.
226 success : `bool`
227 A boolean indicating if all tests have passed.
228 """
229 verifyStats = {}
230 success = True
231 catalogVerify = {}
232 catStats = statisticsDict['CATALOG']
234 catalogVerify['BRIGHT_SLOPE'] = False
235 catalogVerify['NUM_MATCHES'] = False
236 # These values need justification.
237 if catStats['FIT_SUCCESS'] and catStats['BRIGHT_SLOPE'] < -0.5:
238 catalogVerify['BRIGHT_SLOPE'] = True
239 if catStats['NUM_MATCHES'] > 10:
240 catalogVerify['NUM_MATCHES'] = True
242 if catalogVerify['NUM_MATCHES'] and not catalogVerify['BRIGHT_SLOPE']:
243 success = False
244 return {'AMP': verifyStats, 'CATALOG': catalogVerify}, bool(success)
246 def repackStats(self, statisticsDict, dimensions):
247 # docstring inherited
248 rows = {}
249 rowList = []
250 matrixRowList = None
252 rowBase = {
253 "instrument": dimensions["instrument"],
254 "detector": dimensions["detector"],
255 }
257 # CAT results
258 for ampName, stats in statisticsDict["CAT"].items():
259 rows[ampName] = {}
260 rows[ampName].update(rowBase)
261 rows[ampName]["amplifier"] = ampName
262 for key, value in stats.items():
263 rows["ampName"][f"{self.config.stageName}_{key}"] = value
264 rows["ampName"][f"{self.config.stageName}_VERIFY_{key}"] = value
266 # pack final list
267 for ampName, stats in rows.items():
268 rowList.append(stats)
270 return rowList, matrixRowList