Coverage for tests / test_electrostaticBrighterFatter.py: 15%
81 statements
« prev ^ index » next coverage.py v7.13.5, created at 2026-05-05 18:52 +0000
1#
2# LSST Data Management System
3#
4# Copyright 2008-2017 AURA/LSST.
5#
6# This product includes software developed by the
7# LSST Project (http://www.lsst.org/).
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 LSST License Statement and
20# the GNU General Public License along with this program. If not,
21# see:
22# https://www.lsstcorp.org/LegalNotices/
23#
24"""Test cases for lsst.cp.pipe.ElectrostaticBrighterFatterSolveTask.
25"""
27import unittest
28import numpy as np
30import lsst.utils
31import lsst.utils.tests
33import lsst.ip.isr as ipIsr
34import lsst.cp.pipe as cpPipe
35import lsst.afw.cameraGeom as cameraGeom
38class ElectrostaticBrighterFatterSolveTaskTestCase(lsst.utils.tests.TestCase):
39 """ A test case for the electrostatic brighter-fatter solver.
40 """
42 def setUp(self):
43 """ Set up an empirical PTC dataset.
44 """
45 cameraBuilder = cameraGeom.Camera.Builder('fake camera')
46 detectorWrapper = cameraGeom.testUtils.DetectorWrapper(
47 numAmps=16, cameraBuilder=cameraBuilder,
48 )
50 self.detector = detectorWrapper.detector
52 self.camera = cameraBuilder.finish()
54 self.defaultConfig = cpPipe.BrighterFatterKernelSolveConfig()
55 ampNames = [f"amp {i}" for i in range(16)]
57 self.ptc = ipIsr.PhotonTransferCurveDataset(
58 ampNames=ampNames, ptcFitType='FULLCOVARIANCE',
59 covMatrixSide=5,
60 )
62 self.ptc.badAmps = [ampNames[3]] # Randomly set one amp as bad
64 # These comes from the empirical average a matrix from
65 # LSSTCam/calib/DM-51897 for detector R22_S11
66 self.aMatrixMean = np.array([
67 [-3.29036958e-06, 4.16003921e-07, 6.19957768e-08,
68 1.73248403e-08, 1.15894578e-08],
69 [1.75050697e-07, 1.25755064e-07, 4.21671079e-08,
70 1.70534397e-08, 6.97588077e-09],
71 [5.90596810e-08, 3.85312617e-08, 2.25833070e-08,
72 9.42166010e-09, 8.42694139e-09],
73 [1.95655870e-08, 1.74141350e-08, 1.09739771e-08,
74 8.40992104e-09, 3.83679509e-09],
75 [1.10096895e-08, 7.17968238e-09, 7.99943349e-09,
76 2.64435974e-09, 4.21998095e-09],
77 ])
79 self.aMatrixSigma = np.array([
80 [43.61315413, 5.25852337, 3.49826154, 3.25020952,
81 3.70877012],
82 [6.24446168, 2.80054151, 2.36052852, 3.03517635,
83 2.75386129],
84 [4.52060910, 2.87046950, 2.81298062, 3.23891077,
85 3.48839977],
86 [4.78699906, 3.56665023, 2.49834086, 2.93183221,
87 2.58406590],
88 [3.08975509, 2.39008599, 2.76165328, 2.83526835,
89 3.14183071],
90 ]) * 1e-9
92 rawMeans = np.array([
93 1000, 2000, 3000, 4000, 5000,
94 6000, 7000, 8000, 9000, 10000
95 ], dtype=np.float64)
97 rawVars = np.array([
98 708.54575502, 1398.11082581, 2081.09149706,
99 2757.56675448, 3427.61457982, 4091.31195080,
100 4748.73484115, 5399.95822064, 6045.05605499,
101 6684.10130596,
102 ])
104 for i, ampName in enumerate(ampNames):
105 rng = np.random.default_rng(i)
106 self.ptc.aMatrix[ampName] = rng.normal(
107 loc=self.aMatrixMean, scale=self.aMatrixSigma,
108 )
109 self.ptc.gain[ampName] = 1.0
110 self.ptc.rawMeans[ampName] = rawMeans
111 self.ptc.rawVars[ampName] = rawVars
112 self.ptc.finalMeans[ampName] = rawMeans
113 self.ptc.finalVars[ampName] = rawVars
114 self.ptc.expIdMask[ampName] = np.ones_like(rawMeans, dtype=bool)
115 self.ptc.covariances[ampName] = []
116 for mean, variance in zip(rawMeans, rawVars):
117 residual = mean - variance
118 m = np.array([
119 [1, -2e-2, -3e-3, 3e-3, 0e0],
120 [-5e-3, -2e-03, 1e-3, 0e0, 0e0],
121 [-3e-3, -2e-3, 1e-3, 2e-3, 1e-3],
122 [-2e-3, -1e-3, 2e-3, 3e-3, 0e0],
123 [-1e-3, 1e-3, -2e-3, 0e0, 0e0],
124 ]) # 5x5
125 covariance = m * np.full_like(m, residual)
126 covariance[0][0] = variance
128 self.ptc.covariances[ampName].append(covariance)
130 self.ptc.covariancesModel[ampName] = (
131 self.ptc.covariances[ampName]
132 )
133 self.ptc.gain[ampName] = 1.0
134 self.ptc.noise[ampName] = 5.0
135 self.ptc.noiseMatrix[ampName] = np.zeros((5, 5))
136 self.ptc.noiseMatrix[ampName][0][0] = (
137 self.ptc.noise[ampName]**2
138 )
140 self.sequencerMetadata = {
141 "SEQNAME": "a_sequencer",
142 "SEQFILE": "a_sequencer_file",
143 "SEQCKSUM": "hedgehog",
144 }
145 self.ptc.updateMetadata(
146 **self.sequencerMetadata,
147 setCalibInfo=True,
148 )
150 # This is empirically determined from the above parameters.
151 self.expectationN = np.array([
152 [-9.6043345444233e-07, -1.834928531884573e-07,
153 -6.829002404001458e-08, -3.395009268910933e-08,
154 -1.938577698271509e-08],
155 [-1.67723647873235e-07, -1.1170727096154204e-07,
156 -5.4937837964336995e-08, -3.0054681985728955e-08,
157 -1.7905956667356382e-08],
158 [-2.613505487141988e-08, -4.248790474354675e-08,
159 -3.250610285914959e-08, -2.175383668678412e-08,
160 -1.435793965095469e-08],
161 [-7.800476053427953e-09, -1.7070646760313933e-08,
162 -1.7368487910271676e-08, -1.4026417648212689e-08,
163 -1.0402838396107904e-08],
164 [-3.16126749248584e-09, -7.815752671871989e-09,
165 -9.321812906402283e-09, -8.633007891598854e-09,
166 -7.093969018463969e-09]
167 ])
169 self.expectationS = np.array([
170 [-9.604334544423303e-07, 9.6043345444233e-07,
171 1.834928531884573e-07, 6.829002404001458e-08,
172 3.395009268910933e-08],
173 [-1.67723647873235e-07, 1.67723647873235e-07,
174 1.1170727096154204e-07, 5.4937837964336995e-08,
175 3.0054681985728955e-08],
176 [-2.6135054871419876e-08, 2.613505487141988e-08,
177 4.248790474354675e-08, 3.250610285914959e-08,
178 2.175383668678412e-08],
179 [-7.800476053427953e-09, 7.800476053427953e-09,
180 1.7070646760313933e-08, 1.7368487910271676e-08,
181 1.4026417648212689e-08],
182 [-3.16126749248584e-09, 3.16126749248584e-09,
183 7.815752671871989e-09, 9.321812906402283e-09,
184 8.633007891598854e-09]
185 ])
187 self.expectationE = np.array([
188 [-6.749744427777267e-07, -1.8250339728555552e-07,
189 -2.7729149854620917e-08, -8.035079968351e-09,
190 -3.2184346783570578e-09],
191 [-1.6610528501534137e-07, -1.0730697633591896e-07,
192 -4.285764296461245e-08, -1.7319955895167644e-08,
193 -7.91199524741906e-09],
194 [-6.549057579118752e-08, -5.3404306350449054e-08,
195 -3.225504115994437e-08, -1.7412759905918106e-08,
196 -9.374909613003027e-09],
197 [-3.3172336024330024e-08, -2.949619795693681e-08,
198 -2.153854786613826e-08, -1.3986334903747361e-08,
199 -8.642359667681748e-09],
200 [-1.9093161370773104e-08, -1.766815256291286e-08,
201 -1.4227394294305891e-08, -1.03532572194309e-08,
202 -7.0831579167455315e-09]
203 ])
205 self.expectationW = np.array([
206 [-6.749744427777267e-07, -1.8250339728555552e-07,
207 -2.7729149854620917e-08, -8.035079968351e-09,
208 -3.2184346783570578e-09],
209 [6.749744427777267e-07, 1.8250339728555552e-07,
210 2.7729149854620917e-08, 8.035079968351e-09,
211 3.2184346783570578e-09],
212 [1.6610528501534137e-07, 1.0730697633591896e-07,
213 4.285764296461245e-08, 1.7319955895167644e-08,
214 7.91199524741906e-09],
215 [6.549057579118752e-08, 5.3404306350449054e-08,
216 3.225504115994437e-08, 1.7412759905918106e-08,
217 9.374909613003027e-09],
218 [3.3172336024330024e-08, 2.949619795693681e-08,
219 2.153854786613826e-08, 1.3986334903747361e-08,
220 8.642359667681748e-09]
221 ])
223 self.expectationFitParams = {
224 "thickness": 100.0,
225 "pixelsize": 10.0,
226 "zq": 2.90316973,
227 "zsh_minus_zq": 8.2469e-07,
228 "zsh": 2.90317056,
229 "zsv_minus_zq": 1.09756267,
230 "zsv": 4.00073241,
231 "a": 0.00191671,
232 "b": 3.5,
233 "alpha": 0.64506325,
234 "beta": 7.2694e-10,
235 }
237 def test_average(self):
238 """Test "averaged" input aMatrix and aMatrixSigma.
239 """
240 config = cpPipe.ElectrostaticBrighterFatterSolveTask().ConfigClass()
241 config.fitRange = 3
242 task = cpPipe.ElectrostaticBrighterFatterSolveTask(config=config)
244 results = task.run(self.ptc, dummy=['this is a dummy exposure'],
245 camera=self.camera, inputDims={'detector': 1})
247 electroBfDistortionMatrix = results.output
249 self.assertFloatsAlmostEqual(
250 electroBfDistortionMatrix.aMatrix, self.aMatrixMean,
251 atol=3*self.aMatrixSigma,
252 )
253 self.assertFloatsAlmostEqual(
254 electroBfDistortionMatrix.aMatrixSigma, self.aMatrixSigma,
255 atol=3*self.aMatrixSigma/np.sqrt(8)
256 )
258 for key, value in self.sequencerMetadata.items():
259 self.assertEqual(electroBfDistortionMatrix.metadata[key], value)
261 self.assertEqual(electroBfDistortionMatrix.metadata["INSTRUME"], self.camera.getName())
262 self.assertEqual(electroBfDistortionMatrix.metadata["DETECTOR"], self.detector.getId())
263 self.assertEqual(electroBfDistortionMatrix.metadata["DET_NAME"], self.detector.getName())
264 self.assertEqual(electroBfDistortionMatrix.metadata["DET_SER"], self.detector.getSerial())
266 def test_electrostaticSolution(self):
267 """Test the full electrostatic solution.
268 """
269 config = cpPipe.ElectrostaticBrighterFatterSolveTask().ConfigClass()
270 config.fitRange = 5
271 task = cpPipe.ElectrostaticBrighterFatterSolveTask(config=config)
273 results = task.run(self.ptc, dummy=['this is a dummy exposure'],
274 camera=self.camera, inputDims={'detector': 1})
275 electroBfDistortionMatrix = results.output
277 for n in electroBfDistortionMatrix.fitParamNames:
278 self.assertFloatsAlmostEqual(
279 electroBfDistortionMatrix.fitParams[n],
280 self.expectationFitParams[n],
281 atol=1e-3,
282 )
284 d = 1
285 self.assertFloatsAlmostEqual(
286 electroBfDistortionMatrix.aN[:d, :d], self.expectationN[:d, :d],
287 rtol=1e-3,
288 )
289 self.assertFloatsAlmostEqual(
290 electroBfDistortionMatrix.aS[:d, :d], self.expectationS[:d, :d],
291 rtol=1e-3,
292 )
293 self.assertFloatsAlmostEqual(
294 electroBfDistortionMatrix.aE[:d, :d], self.expectationE[:d, :d],
295 rtol=1e-3,
296 )
297 self.assertFloatsAlmostEqual(
298 electroBfDistortionMatrix.aW[:d, :d], self.expectationW[:d, :d],
299 rtol=1e-3,
300 )
303def setup_module(module):
304 lsst.utils.tests.init()
307if __name__ == "__main__": 307 ↛ 308line 307 didn't jump to line 308 because the condition on line 307 was never true
308 lsst.utils.tests.init()
309 unittest.main()