Coverage for tests / test_electrostaticBrighterFatter.py: 15%

82 statements  

« prev     ^ index     » next       coverage.py v7.13.5, created at 2026-04-28 09:20 +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""" 

26 

27import unittest 

28import numpy as np 

29 

30import lsst.utils 

31import lsst.utils.tests 

32 

33import lsst.ip.isr as ipIsr 

34import lsst.cp.pipe as cpPipe 

35import lsst.afw.cameraGeom as cameraGeom 

36 

37 

38class ElectrostaticBrighterFatterSolveTaskTestCase(lsst.utils.tests.TestCase): 

39 """ A test case for the electrostatic brighter-fatter solver. 

40 """ 

41 

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 ) 

49 

50 self.detector = detectorWrapper.detector 

51 

52 self.camera = cameraBuilder.finish() 

53 

54 self.defaultConfig = cpPipe.BrighterFatterKernelSolveConfig() 

55 ampNames = [f"amp {i}" for i in range(16)] 

56 

57 self.ptc = ipIsr.PhotonTransferCurveDataset( 

58 ampNames=ampNames, ptcFitType='FULLCOVARIANCE', 

59 covMatrixSide=5, 

60 ) 

61 

62 self.ptc.badAmps = [ampNames[3]] # Randomly set one amp as bad 

63 

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 ]) 

78 

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 

91 

92 rawMeans = np.array([ 

93 1000, 2000, 3000, 4000, 5000, 

94 6000, 7000, 8000, 9000, 10000 

95 ], dtype=np.float64) 

96 

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 ]) 

103 

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 

127 

128 self.ptc.covariances[ampName].append(covariance) 

129 

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 ) 

139 

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 ) 

149 

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 ]) 

168 

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 ]) 

186 

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 ]) 

204 

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 ]) 

222 

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 } 

236 

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) 

243 

244 results = task.run(self.ptc, dummy=['this is a dummy exposure'], 

245 camera=self.camera, inputDims={'detector': 1}) 

246 

247 electroBfDistortionMatrix = results.output 

248 

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 ) 

257 

258 for key, value in self.sequencerMetadata.items(): 

259 self.assertEqual(electroBfDistortionMatrix.metadata[key], value) 

260 

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()) 

265 

266 def test_electrostaticSolution(self): 

267 """Test the full electrostatic solution. 

268 """ 

269 config = cpPipe.ElectrostaticBrighterFatterSolveTask().ConfigClass() 

270 config.fitRange = 5 

271 config.doFilterCorrection = False 

272 task = cpPipe.ElectrostaticBrighterFatterSolveTask(config=config) 

273 

274 results = task.run(self.ptc, dummy=["This is a dummy exposure"], 

275 camera=self.camera, inputDims={'detector': 1}) 

276 electroBfDistortionMatrix = results.output 

277 

278 for n in electroBfDistortionMatrix.fitParamNames: 

279 self.assertFloatsAlmostEqual( 

280 electroBfDistortionMatrix.fitParams[n], 

281 self.expectationFitParams[n], 

282 atol=1e-3, 

283 ) 

284 

285 d = 1 

286 self.assertFloatsAlmostEqual( 

287 electroBfDistortionMatrix.aN[:d, :d], self.expectationN[:d, :d], 

288 rtol=1e-3, 

289 ) 

290 self.assertFloatsAlmostEqual( 

291 electroBfDistortionMatrix.aS[:d, :d], self.expectationS[:d, :d], 

292 rtol=1e-3, 

293 ) 

294 self.assertFloatsAlmostEqual( 

295 electroBfDistortionMatrix.aE[:d, :d], self.expectationE[:d, :d], 

296 rtol=1e-3, 

297 ) 

298 self.assertFloatsAlmostEqual( 

299 electroBfDistortionMatrix.aW[:d, :d], self.expectationW[:d, :d], 

300 rtol=1e-3, 

301 ) 

302 

303 

304def setup_module(module): 

305 lsst.utils.tests.init() 

306 

307 

308if __name__ == "__main__": 308 ↛ 309line 308 didn't jump to line 309 because the condition on line 308 was never true

309 lsst.utils.tests.init() 

310 unittest.main()