Coverage for tests / test_statistics.py: 14%

361 statements  

« prev     ^ index     » next       coverage.py v7.13.5, created at 2026-04-17 02:05 -0700

1# This file is part of afw. 

2# 

3# Developed for the LSST Data Management System. 

4# This product includes software developed by the LSST Project 

5# (https://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 <https://www.gnu.org/licenses/>. 

21 

22""" 

23Tests for Statistics 

24 

25Run with: 

26 python test_statistics.py 

27or 

28 pytest test_statistics.py 

29""" 

30 

31import math 

32import os 

33import unittest 

34 

35import numpy as np 

36 

37import lsst.utils.tests 

38import lsst.pex.exceptions 

39import lsst.geom 

40import lsst.afw.image as afwImage 

41import lsst.afw.math as afwMath 

42import lsst.afw.display as afwDisplay 

43 

44afwDisplay.setDefaultMaskTransparency(75) 

45 

46try: 

47 afwdataDir = lsst.utils.getPackageDir("afwdata") 

48except LookupError: 

49 afwdataDir = None 

50 

51try: 

52 type(display) 

53except NameError: 

54 display = False 

55 

56 

57class StatisticsTestCase(lsst.utils.tests.TestCase): 

58 """A test case for Statistics""" 

59 

60 clippedVariance3 = 0.9733369 # variance of an N(0, 1) Gaussian clipped at 3 sigma 

61 

62 def setUp(self): 

63 w, h = 900, 1500 

64 

65 mean = 10.5 # requested mean 

66 std = 1.0 # and standard deviation 

67 

68 self.images = [] 

69 # ImageI 

70 np.random.seed(666) 

71 isInt = True 

72 image = afwImage.ImageI(lsst.geom.ExtentI(w, h)) 

73 image.array[:] = np.floor(np.random.normal(mean, std, (h, w)) + 0.5).astype(int) 

74 

75 # Note that the mean/median/std may not be quite equal to the requested values 

76 self.images.append((image, isInt, np.mean(image.array), np.mean(image.array), np.std(image.array))) 

77 

78 # ImageF 

79 np.random.seed(666) 

80 isInt = False 

81 image = afwImage.ImageF(lsst.geom.ExtentI(w, h)) 

82 image.array[:] = np.random.normal(mean, std, (h, w)) 

83 

84 # Note that the mean/median/std may not be quite equal to the requested values 

85 self.images.append((image, isInt, np.mean(image.array), np.median(image.array), np.std(image.array))) 

86 

87 @staticmethod 

88 def delta(what, isInt): 

89 # Return a tolerance for a test 

90 if what == "mean": 

91 return 4e-6 

92 elif what == "meanclip": 

93 return 4e-5 

94 elif what == "median": 

95 return 0.00022 if isInt else 0.00000075 

96 

97 def tearDown(self): 

98 del self.images 

99 

100 def testDefaultGet(self): 

101 """Test that we can get a single statistic without specifying it""" 

102 for image, isInt, mean, median, std in self.images: 

103 stats = afwMath.makeStatistics(image, afwMath.MEDIAN) 

104 

105 self.assertEqual(stats.getValue(), stats.getValue(afwMath.MEDIAN)) 

106 self.assertEqual(stats.getResult()[0], stats.getResult(afwMath.MEDIAN)[0]) 

107 # 

108 stats = afwMath.makeStatistics(image, afwMath.MEDIAN | afwMath.ERRORS) 

109 

110 self.assertEqual(stats.getValue(), stats.getValue(afwMath.MEDIAN)) 

111 self.assertEqual(stats.getResult(), stats.getResult(afwMath.MEDIAN)) 

112 self.assertEqual(stats.getError(), stats.getError(afwMath.MEDIAN)) 

113 

114 def tst(): 

115 stats.getValue() 

116 

117 stats = afwMath.makeStatistics(image, afwMath.MEDIAN | afwMath.MEAN) 

118 self.assertRaises(lsst.pex.exceptions.InvalidParameterError, tst) 

119 

120 def testStats1(self): 

121 for image, isInt, mean, median, std in self.images: 

122 stats = afwMath.makeStatistics(image, afwMath.NPOINT | afwMath.STDEV | afwMath.MEAN | afwMath.SUM) 

123 

124 self.assertEqual(stats.getValue(afwMath.NPOINT), image.getWidth() * image.getHeight()) 

125 self.assertEqual( 

126 stats.getValue(afwMath.NPOINT) * stats.getValue(afwMath.MEAN), stats.getValue(afwMath.SUM) 

127 ) 

128 

129 self.assertAlmostEqual(stats.getValue(afwMath.MEAN), mean, delta=self.delta("mean", isInt)) 

130 # didn't ask for error, so it's a NaN 

131 self.assertTrue(np.isnan(stats.getError(afwMath.MEAN))) 

132 self.assertAlmostEqual(stats.getValue(afwMath.STDEV), std, delta=0.000008) 

133 

134 def testStats2(self): 

135 for image, isInt, mean, median, std in self.images: 

136 stats = afwMath.makeStatistics(image, afwMath.STDEV | afwMath.MEAN | afwMath.ERRORS) 

137 meanRes = stats.getResult(afwMath.MEAN) 

138 sd = stats.getValue(afwMath.STDEV) 

139 

140 self.assertAlmostEqual(meanRes[0], mean, delta=self.delta("mean", isInt)) 

141 self.assertAlmostEqual(meanRes[1], sd / math.sqrt(image.getWidth() * image.getHeight())) 

142 

143 def testStats3(self): 

144 for image, isInt, mean, median, std in self.images: 

145 stats = afwMath.makeStatistics(image, afwMath.NPOINT) 

146 

147 def getMean(): 

148 stats.getValue(afwMath.MEAN) 

149 

150 self.assertRaises(lsst.pex.exceptions.InvalidParameterError, getMean) 

151 

152 def testStatsZebra(self): 

153 """Add 1 to every other row""" 

154 for image, isInt, mean, median, std in self.images: 

155 image2 = image.clone() 

156 # 

157 # Add 1 to every other row, so the variance is increased by 1/4 

158 # 

159 self.assertEqual(image2.getHeight() % 2, 0) 

160 width = image2.getWidth() 

161 for y in range(1, image2.getHeight(), 2): 

162 sim = image2[lsst.geom.Box2I(lsst.geom.Point2I(0, y), lsst.geom.Extent2I(width, 1))] 

163 sim += 1 

164 

165 if display: 

166 afwDisplay.Display(frame=0).mtv(image, "Image 1") 

167 afwDisplay.Display(frame=1).mtv(image2, "Image 2 (var inc by 1/4)") 

168 

169 stats = afwMath.makeStatistics( 

170 image2, afwMath.NPOINT | afwMath.STDEV | afwMath.MEAN | afwMath.ERRORS 

171 ) 

172 meanRes = stats.getResult(afwMath.MEAN) 

173 n = stats.getValue(afwMath.NPOINT) 

174 sd = stats.getValue(afwMath.STDEV) 

175 

176 self.assertAlmostEqual(meanRes[0], mean + 0.5, delta=self.delta("mean", isInt)) 

177 self.assertAlmostEqual( 

178 sd, np.hypot(std, 1 / math.sqrt(4.0) * math.sqrt(n / (n - 1))), delta=0.00011 

179 ) 

180 self.assertAlmostEqual(meanRes[1], sd / math.sqrt(image2.getWidth() * image2.getHeight()), 10) 

181 

182 meanSquare = afwMath.makeStatistics(image2, afwMath.MEANSQUARE).getValue() 

183 self.assertAlmostEqual( 

184 meanSquare, 0.5 * (mean**2 + (mean + 1) ** 2) + std**2, delta=0.00025 if isInt else 0.0005 

185 ) 

186 

187 def testStatsStdevclip(self): 

188 """Test STDEVCLIP; cf. #611""" 

189 for image, isInt, mean, median, std in self.images: 

190 image2 = image.clone() 

191 

192 stats = afwMath.makeStatistics(image2, afwMath.STDEVCLIP | afwMath.NPOINT | afwMath.SUM) 

193 self.assertAlmostEqual( 

194 stats.getValue(afwMath.STDEVCLIP), math.sqrt(self.clippedVariance3) * std, delta=0.0015 

195 ) 

196 # 

197 # Check we get the correct sum even when clipping 

198 # 

199 self.assertEqual( 

200 stats.getValue(afwMath.NPOINT) * afwMath.makeStatistics(image2, afwMath.MEAN).getValue(), 

201 stats.getValue(afwMath.SUM), 

202 ) 

203 

204 def testMedian(self): 

205 """Test the median code""" 

206 for image, isInt, mean, median, std in self.images: 

207 med = afwMath.makeStatistics(image, afwMath.MEDIAN).getValue() 

208 self.assertAlmostEqual(med, median, delta=self.delta("median", isInt)) 

209 

210 values = [1.0, 2.0, 3.0, 2.0] 

211 self.assertEqual(afwMath.makeStatistics(values, afwMath.MEDIAN).getValue(), 2.0) 

212 

213 def testIqrange(self): 

214 """Test the inter-quartile range""" 

215 for image, isInt, mean, median, std in self.images: 

216 iqr = afwMath.makeStatistics(image, afwMath.IQRANGE).getValue() 

217 # pretty loose constraint for isInt; probably because the distribution 

218 # isn't very Gaussian with the added rounding to integer values 

219 self.assertAlmostEqual(iqr, std / 0.741301109252802, delta=0.063 if isInt else 0.00011) 

220 

221 def testMeanClip(self): 

222 """Test the clipped mean""" 

223 

224 sctrl = afwMath.StatisticsControl() 

225 sctrl.setNumSigmaClip(6) 

226 

227 for image, isInt, mean, median, std in self.images: 

228 stats = afwMath.makeStatistics(image, afwMath.MEANCLIP | afwMath.NCLIPPED, sctrl) 

229 self.assertAlmostEqual(stats.getValue(afwMath.MEANCLIP), mean, delta=self.delta("mean", isInt)) 

230 self.assertEqual(stats.getValue(afwMath.NCLIPPED), 0) 

231 

232 def testVarianceClip(self): 

233 """Test the 3-sigma clipped standard deviation and variance""" 

234 for image, isInt, mean, median, std in self.images: 

235 delta = 0.0006 if isInt else 0.0014 

236 stdevClip = afwMath.makeStatistics(image, afwMath.STDEVCLIP).getValue() 

237 self.assertAlmostEqual(stdevClip, math.sqrt(self.clippedVariance3) * std, delta=delta) 

238 

239 varianceClip = afwMath.makeStatistics(image, afwMath.VARIANCECLIP).getValue() 

240 self.assertAlmostEqual(varianceClip, self.clippedVariance3 * std**2, delta=2 * delta) 

241 

242 def _testBadValue(self, badVal): 

243 """Test that we can handle an instance of `badVal` in the data correctly 

244 

245 Note that we only test ImageF here (as ImageI can't contain a NaN) 

246 """ 

247 mean = self.images[0][1] 

248 x, y = 10, 10 

249 for useImage in [True, False]: 

250 if useImage: 

251 image = afwImage.ImageF(100, 100) 

252 image.set(mean) 

253 image[x, y] = badVal 

254 else: 

255 image = afwImage.MaskedImageF(100, 100) 

256 image.set(mean, 0x0, 1.0) 

257 image[x, y] = (badVal, 0x0, 1.0) 

258 

259 self.assertEqual(afwMath.makeStatistics(image, afwMath.MAX).getValue(), mean) 

260 self.assertEqual(afwMath.makeStatistics(image, afwMath.MEAN).getValue(), mean) 

261 

262 sctrl = afwMath.StatisticsControl() 

263 

264 sctrl.setNanSafe(False) 

265 self.assertFalse(np.isfinite(afwMath.makeStatistics(image, afwMath.MAX, sctrl).getValue())) 

266 self.assertFalse(np.isfinite(afwMath.makeStatistics(image, afwMath.MEAN, sctrl).getValue())) 

267 

268 def testMaxWithNan(self): 

269 """Test that we can handle NaNs correctly""" 

270 self._testBadValue(np.nan) 

271 

272 def testMaxWithInf(self): 

273 """Test that we can handle infinities correctly""" 

274 self._testBadValue(np.inf) 

275 

276 @unittest.skipIf(afwdataDir is None, "afwdata not setup") 

277 def testSampleImageStats(self): 

278 """Compare our results to known values in test data""" 

279 

280 imgfiles = [] 

281 imgfiles.append("v1_i1_g_m400_s20_f.fits") 

282 imgfiles.append("v1_i1_g_m400_s20_u16.fits") 

283 imgfiles.append("v1_i2_g_m400_s20_f.fits") 

284 imgfiles.append("v1_i2_g_m400_s20_u16.fits") 

285 imgfiles.append("v2_i1_p_m9_f.fits") 

286 imgfiles.append("v2_i1_p_m9_u16.fits") 

287 imgfiles.append("v2_i2_p_m9_f.fits") 

288 imgfiles.append("v2_i2_p_m9_u16.fits") 

289 

290 afwdataDir = os.getenv("AFWDATA_DIR") 

291 

292 for imgfile in imgfiles: 

293 imgPath = os.path.join(afwdataDir, "Statistics", imgfile) 

294 

295 # get the image and header 

296 dimg = afwImage.DecoratedImageF(imgPath) 

297 fitsHdr = dimg.getMetadata() 

298 

299 # get the true values of the mean and stdev 

300 trueMean = fitsHdr.getAsDouble("MEANCOMP") 

301 trueStdev = fitsHdr.getAsDouble("SIGCOMP") 

302 

303 # measure the mean and stdev with the Statistics class 

304 img = dimg.getImage() 

305 statobj = afwMath.makeStatistics(img, afwMath.MEAN | afwMath.STDEV) 

306 mean = statobj.getValue(afwMath.MEAN) 

307 stdev = statobj.getValue(afwMath.STDEV) 

308 

309 # print trueMean, mean, trueStdev, stdev 

310 self.assertAlmostEqual(mean, trueMean, 8) 

311 self.assertAlmostEqual(stdev, trueStdev, 8) 

312 

313 def testStatisticsRamp(self): 

314 """Tests Statistics on a 'ramp' (image with constant gradient)""" 

315 

316 nx = 101 

317 ny = 64 

318 img = afwImage.ImageF(lsst.geom.Extent2I(nx, ny)) 

319 

320 z0 = 10.0 

321 dzdx = 1.0 

322 mean = z0 + (nx // 2) * dzdx 

323 stdev = 0.0 

324 for y in range(ny): 

325 for x in range(nx): 

326 z = z0 + dzdx * x 

327 img[x, y] = z 

328 stdev += (z - mean) * (z - mean) 

329 

330 stdev = math.sqrt(stdev / (nx * ny - 1)) 

331 

332 stats = afwMath.makeStatistics(img, afwMath.NPOINT | afwMath.STDEV | afwMath.MEAN) 

333 testmean = stats.getValue(afwMath.MEAN) 

334 teststdev = stats.getValue(afwMath.STDEV) 

335 

336 self.assertEqual(stats.getValue(afwMath.NPOINT), nx * ny) 

337 self.assertEqual(testmean, mean) 

338 self.assertAlmostEqual(teststdev, stdev) 

339 

340 stats = afwMath.makeStatistics(img, afwMath.STDEV | afwMath.MEAN | afwMath.ERRORS) 

341 mean, meanErr = stats.getResult(afwMath.MEAN) 

342 sd = stats.getValue(afwMath.STDEV) 

343 

344 self.assertEqual(mean, img[nx // 2, ny // 2]) 

345 self.assertEqual(meanErr, sd / math.sqrt(img.getWidth() * img.getHeight())) 

346 

347 # =============================================================================== 

348 # sjb code for percentiles and clipped stats 

349 

350 stats = afwMath.makeStatistics(img, afwMath.MEDIAN) 

351 self.assertEqual(z0 + dzdx * (nx - 1) / 2.0, stats.getValue(afwMath.MEDIAN)) 

352 

353 stats = afwMath.makeStatistics(img, afwMath.IQRANGE) 

354 self.assertEqual(dzdx * (nx - 1) / 2.0, stats.getValue(afwMath.IQRANGE)) 

355 

356 stats = afwMath.makeStatistics(img, afwMath.MEANCLIP) 

357 self.assertEqual(z0 + dzdx * (nx - 1) / 2.0, stats.getValue(afwMath.MEANCLIP)) 

358 

359 def testMask(self): 

360 mask = afwImage.Mask(lsst.geom.Extent2I(10, 10)) 

361 mask.set(0x0) 

362 

363 mask[1, 1] = 0x10 

364 mask[3, 1] = 0x08 

365 mask[5, 4] = 0x08 

366 mask[4, 5] = 0x02 

367 

368 stats = afwMath.makeStatistics(mask, afwMath.SUM | afwMath.NPOINT) 

369 self.assertEqual(mask.getWidth() * mask.getHeight(), stats.getValue(afwMath.NPOINT)) 

370 self.assertEqual(0x1A, stats.getValue(afwMath.SUM)) 

371 

372 def tst(): 

373 afwMath.makeStatistics(mask, afwMath.MEAN) 

374 

375 self.assertRaises(lsst.pex.exceptions.InvalidParameterError, tst) 

376 

377 def testTicket1025(self): 

378 """ 

379 Ticket #1025 reported that the Statistics median was getting '3' as the median of [1,2,3,2] 

380 it was caused by an off-by-one error in the implementation 

381 """ 

382 

383 # check the exact example in the ticket 

384 values = [1.0, 2.0, 3.0, 2.0] 

385 self.assertEqual(afwMath.makeStatistics(values, afwMath.MEDIAN).getValue(), 2) 

386 self.assertEqual(afwMath.makeStatistics(sorted(values), afwMath.MEDIAN).getValue(), 2) 

387 

388 # check some other possible ways it could show up 

389 values = list(range(10)) 

390 self.assertEqual(afwMath.makeStatistics(values, afwMath.MEDIAN).getValue(), 4.5) 

391 values = list(range(11)) 

392 self.assertEqual(afwMath.makeStatistics(values, afwMath.MEDIAN).getValue(), 5.0) 

393 

394 def testTicket1123(self): 

395 """ 

396 Ticket #1123 reported that the Statistics stack routine throws an exception 

397 when all pixels in a stack are masked. Returning a NaN pixel in the stack is preferred 

398 """ 

399 

400 mean = self.images[0][1] 

401 

402 ctrl = afwMath.StatisticsControl() 

403 ctrl.setAndMask(~0x0) 

404 

405 mimg = afwImage.MaskedImageF(lsst.geom.Extent2I(10, 10)) 

406 mimg.set([mean, 0x1, mean]) 

407 

408 # test the case with no valid pixels ... both mean and stdev should be 

409 # nan 

410 stat = afwMath.makeStatistics(mimg, afwMath.MEAN | afwMath.STDEV, ctrl) 

411 mean = stat.getValue(afwMath.MEAN) 

412 stdev = stat.getValue(afwMath.STDEV) 

413 self.assertNotEqual(mean, mean) # NaN does not equal itself 

414 self.assertNotEqual(stdev, stdev) # NaN does not equal itself 

415 

416 # test the case with one valid pixel ... mean is ok, but stdev should 

417 # still be nan 

418 mimg.getMask()[1, 1] = 0x0 

419 stat = afwMath.makeStatistics(mimg, afwMath.MEAN | afwMath.STDEV, ctrl) 

420 mean = stat.getValue(afwMath.MEAN) 

421 stdev = stat.getValue(afwMath.STDEV) 

422 self.assertEqual(mean, mean) 

423 self.assertNotEqual(stdev, stdev) # NaN does not equal itself 

424 

425 # test the case with two valid pixels ... both mean and stdev are ok 

426 mimg.getMask()[1, 2] = 0x0 

427 stat = afwMath.makeStatistics(mimg, afwMath.MEAN | afwMath.STDEV, ctrl) 

428 mean = stat.getValue(afwMath.MEAN) 

429 stdev = stat.getValue(afwMath.STDEV) 

430 self.assertEqual(mean, mean) 

431 self.assertEqual(stdev, 0.0) 

432 

433 def testTicket1125(self): 

434 """Ticket 1125 reported that the clipped routines were aborting when called with no valid pixels.""" 

435 

436 mean = self.images[0][1] 

437 

438 mimg = afwImage.MaskedImageF(lsst.geom.Extent2I(10, 10)) 

439 mimg.set([mean, 0x1, mean]) 

440 

441 ctrl = afwMath.StatisticsControl() 

442 ctrl.setAndMask(~0x0) 

443 

444 # test the case with no valid pixels ... try MEANCLIP and STDEVCLIP 

445 stat = afwMath.makeStatistics(mimg, afwMath.MEANCLIP | afwMath.STDEVCLIP, ctrl) 

446 mean = stat.getValue(afwMath.MEANCLIP) 

447 stdev = stat.getValue(afwMath.STDEVCLIP) 

448 self.assertNotEqual(mean, mean) # NaN does not equal itself 

449 self.assertNotEqual(stdev, stdev) # NaN does not equal itself 

450 

451 def testWeightedSum(self): 

452 ctrl = afwMath.StatisticsControl() 

453 mi = afwImage.MaskedImageF(lsst.geom.Extent2I(10, 10)) 

454 mi.getImage().set(1.0) 

455 mi.getVariance().set(0.1) 

456 

457 stats = afwMath.makeStatistics(mi, afwMath.SUM, ctrl) 

458 self.assertEqual(stats.getValue(afwMath.SUM), 100.0) 

459 

460 ctrl.setWeighted(True) 

461 weighted = afwMath.makeStatistics(mi, afwMath.SUM, ctrl) 

462 # precision at "4 places" as images are floats 

463 # ... variance = 0.1 is stored as 0.100000001 

464 self.assertAlmostEqual(weighted.getValue(afwMath.SUM), 1000.0, 4) 

465 

466 def testWeightedSum2(self): 

467 """Test using a weight image separate from the variance plane""" 

468 weight, mean = 0.1, 1.0 

469 

470 ctrl = afwMath.StatisticsControl() 

471 mi = afwImage.MaskedImageF(lsst.geom.Extent2I(10, 10)) 

472 npix = 10 * 10 

473 mi.getImage().set(mean) 

474 mi.getVariance().set(np.nan) 

475 

476 weights = afwImage.ImageF(mi.getDimensions()) 

477 weights.set(weight) 

478 

479 stats = afwMath.makeStatistics(mi, afwMath.SUM, ctrl) 

480 self.assertEqual(stats.getValue(afwMath.SUM), mean * npix) 

481 

482 weighted = afwMath.makeStatistics(mi, weights, afwMath.SUM, ctrl) 

483 # precision at "4 places" as images are floats 

484 # ... variance = 0.1 is stored as 0.100000001 

485 self.assertAlmostEqual(weighted.getValue(afwMath.SUM), mean * npix * weight, 4) 

486 

487 def testErrorsFromVariance(self): 

488 """Test that we can estimate the errors from the incoming variances""" 

489 weight, mean, variance = 0.1, 1.0, 10.0 

490 

491 mi = afwImage.MaskedImageF(lsst.geom.Extent2I(10, 10)) 

492 npix = 10 * 10 

493 mi.getImage().set(mean) 

494 mi.getVariance().set(variance) 

495 

496 weights = afwImage.ImageF(mi.getDimensions()) 

497 weights.set(weight) 

498 

499 ctrl = afwMath.StatisticsControl() 

500 ctrl.setCalcErrorFromInputVariance(True) 

501 weighted = afwMath.makeStatistics( 

502 mi, weights, afwMath.MEAN | afwMath.MEANCLIP | afwMath.SUM | afwMath.ERRORS, ctrl 

503 ) 

504 

505 self.assertAlmostEqual(weighted.getValue(afwMath.SUM) / (npix * mean * weight), 1) 

506 self.assertAlmostEqual(weighted.getValue(afwMath.MEAN), mean) 

507 self.assertAlmostEqual(weighted.getError(afwMath.MEAN) ** 2, variance / npix) 

508 self.assertAlmostEqual(weighted.getError(afwMath.MEANCLIP) ** 2, variance / npix) 

509 

510 # calcErrorMosaicMode and calcErrorFromInputVariance cannot both be selected 

511 ctrl.setCalcErrorMosaicMode(True) 

512 with self.assertRaises(lsst.pex.exceptions.InvalidParameterError): 

513 afwMath.makeStatistics(mi, weights, afwMath.MEAN, ctrl) 

514 

515 # Test only mosaic mode wherein output variance is the mean of the input variance 

516 ctrl.setCalcErrorFromInputVariance(False) 

517 statsMosaic = afwMath.makeStatistics( 

518 mi, weights, afwMath.MEAN | afwMath.MEANCLIP | afwMath.SUM | afwMath.ERRORS, ctrl 

519 ) 

520 

521 self.assertAlmostEqual(statsMosaic.getValue(afwMath.SUM) / (npix * mean * weight), 1) 

522 self.assertAlmostEqual(statsMosaic.getValue(afwMath.MEAN), mean) 

523 self.assertAlmostEqual(statsMosaic.getError(afwMath.MEAN) ** 2, variance) 

524 self.assertAlmostEqual(statsMosaic.getError(afwMath.MEANCLIP) ** 2, variance) 

525 

526 def testMeanClipSingleValue(self): 

527 """Verify that the clipped mean doesn't not return NaN for a single value.""" 

528 sctrl = afwMath.StatisticsControl() 

529 sctrl.setNumSigmaClip(6) 

530 

531 for image, isInt, mean, median, std in self.images: 

532 stats = afwMath.makeStatistics(image, afwMath.MEANCLIP | afwMath.NCLIPPED, sctrl) 

533 self.assertAlmostEqual( 

534 stats.getValue(afwMath.MEANCLIP), mean, delta=self.delta("meanclip", isInt) 

535 ) 

536 self.assertEqual(stats.getValue(afwMath.NCLIPPED), 0) 

537 

538 # this bug was caused by the iterative nature of the MEANCLIP. 

539 # With only one point, the sample variance returns NaN to avoid a divide by zero error 

540 # Thus, on the second iteration, the clip width (based on _variance) is NaN and corrupts 

541 # all further calculations. 

542 img = afwImage.ImageF(lsst.geom.Extent2I(1, 1)) 

543 img.set(0) 

544 stats = afwMath.makeStatistics(img, afwMath.MEANCLIP | afwMath.NCLIPPED) 

545 self.assertEqual(stats.getValue(afwMath.MEANCLIP), 0) 

546 self.assertEqual(stats.getValue(afwMath.NCLIPPED), 0) 

547 

548 def testMismatch(self): 

549 """Test that we get an exception when there's a size mismatch""" 

550 scale = 5 

551 for image, isInt, mean, median, std in self.images: 

552 dims = image.getDimensions() 

553 mask = afwImage.Mask(dims * scale) 

554 mask.set(0xFF) 

555 ctrl = afwMath.StatisticsControl() 

556 ctrl.setAndMask(0xFF) 

557 # If it didn't raise, this would result in a NaN (the image data is 

558 # completely masked). 

559 self.assertRaises( 

560 lsst.pex.exceptions.InvalidParameterError, 

561 afwMath.makeStatistics, 

562 image, 

563 mask, 

564 afwMath.MEDIAN, 

565 ctrl, 

566 ) 

567 subMask = afwImage.Mask(mask, lsst.geom.Box2I(lsst.geom.Point2I(dims * (scale - 1)), dims)) 

568 subMask.set(0) 

569 # Using subMask is successful. 

570 self.assertAlmostEqual( 

571 afwMath.makeStatistics(image, subMask, afwMath.MEDIAN, ctrl).getValue(), 

572 median, 

573 delta=self.delta("median", isInt), 

574 ) 

575 

576 def testClipping(self): 

577 """Test that clipping statistics work 

578 

579 Insert a single bad pixel; it should be clipped. 

580 """ 

581 sctrl = afwMath.StatisticsControl() 

582 sctrl.setNumSigmaClip(10) 

583 

584 for image, isInt, mean, median, std in self.images: 

585 nval = 1000 * mean 

586 if isInt: 

587 nval = int(nval) 

588 image[0, 0] = nval 

589 

590 stats = afwMath.makeStatistics(image, afwMath.MEANCLIP | afwMath.NCLIPPED | afwMath.NPOINT, sctrl) 

591 self.assertAlmostEqual( 

592 stats.getValue(afwMath.MEANCLIP), mean, delta=self.delta("meanclip", isInt) 

593 ) 

594 self.assertEqual(stats.getValue(afwMath.NCLIPPED), 1) 

595 self.assertEqual(stats.getValue(afwMath.NPOINT), image.getBBox().getArea()) 

596 

597 def testNMasked(self): 

598 """Test that NMASKED works""" 

599 maskVal = 0xBE 

600 ctrl = afwMath.StatisticsControl() 

601 ctrl.setAndMask(maskVal) 

602 for image, isInt, mean, median, std in self.images: 

603 mask = afwImage.Mask(image.getBBox()) 

604 mask.set(0) 

605 self.assertEqual(afwMath.makeStatistics(image, mask, afwMath.NMASKED, ctrl).getValue(), 0) 

606 mask[1, 1] = maskVal 

607 self.assertEqual(afwMath.makeStatistics(image, mask, afwMath.NMASKED, ctrl).getValue(), 1) 

608 

609 

610class TestMemory(lsst.utils.tests.MemoryTestCase): 

611 pass 

612 

613 

614def setup_module(module): 

615 lsst.utils.tests.init() 

616 

617 

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

619 lsst.utils.tests.init() 

620 unittest.main()