Coverage for tests / test_sourceTable.py: 10%

553 statements  

« prev     ^ index     » next       coverage.py v7.13.1, created at 2026-01-13 01:27 -0800

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 

22import os 

23import unittest 

24import tempfile 

25import pickle 

26import math 

27 

28import numpy as np 

29 

30import lsst.utils.tests 

31import lsst.pex.exceptions 

32import lsst.geom 

33import lsst.afw.table 

34import lsst.afw.geom 

35import lsst.afw.image 

36import lsst.afw.detection 

37 

38testPath = os.path.abspath(os.path.dirname(__file__)) 

39 

40 

41def makeArray(size, dtype): 

42 return np.array(np.random.randn(*size), dtype=dtype) 

43 

44 

45def makeCov(size, dtype): 

46 m = np.array(np.random.randn(size, size), dtype=dtype) 

47 return np.dot(m, m.transpose()) 

48 

49 

50def makeWcs(): 

51 crval = lsst.geom.SpherePoint(1.606631*lsst.geom.degrees, 

52 5.090329*lsst.geom.degrees) 

53 crpix = lsst.geom.Point2D(2036.0, 2000.0) 

54 cdMatrix = np.array([5.399452e-5, -1.30770e-5, 1.30770e-5, 5.399452e-5]) 

55 cdMatrix.shape = (2, 2) 

56 return lsst.afw.geom.makeSkyWcs(crpix=crpix, crval=crval, cdMatrix=cdMatrix) 

57 

58 

59class SourceTableTestCase(lsst.utils.tests.TestCase): 

60 

61 def fillRecord(self, record): 

62 record.set(self.instFluxKey, np.random.randn()) 

63 record.set(self.instFluxErrKey, np.random.randn()) 

64 record.set(self.centroidKey.getX(), np.random.randn()) 

65 record.set(self.centroidKey.getY(), np.random.randn()) 

66 record.set(self.xErrKey, np.random.randn()) 

67 record.set(self.yErrKey, np.random.randn()) 

68 record.set(self.shapeKey.getIxx(), np.random.randn()) 

69 record.set(self.shapeKey.getIyy(), np.random.randn()) 

70 record.set(self.shapeKey.getIxy(), np.random.randn()) 

71 record.set(self.xxErrKey, np.random.randn()) 

72 record.set(self.yyErrKey, np.random.randn()) 

73 record.set(self.xyErrKey, np.random.randn()) 

74 record.set(self.psfShapeKey.getIxx(), np.random.randn()) 

75 record.set(self.psfShapeKey.getIyy(), np.random.randn()) 

76 record.set(self.psfShapeKey.getIxy(), np.random.randn()) 

77 record.set(self.fluxFlagKey, np.random.randn() > 0) 

78 record.set(self.centroidFlagKey, np.random.randn() > 0) 

79 record.set(self.shapeFlagKey, np.random.randn() > 0) 

80 

81 def setUp(self): 

82 np.random.seed(1) 

83 self.schema = lsst.afw.table.SourceTable.makeMinimalSchema() 

84 self.instFluxKey = self.schema.addField("a_instFlux", type="D") 

85 self.instFluxErrKey = self.schema.addField("a_instFluxErr", type="D") 

86 self.fluxFlagKey = self.schema.addField("a_flag", type="Flag") 

87 

88 # the meas field is added using a functor key, but the error is added 

89 # as scalars, as we lack a ResultKey functor as exists in meas_base 

90 self.centroidKey = lsst.afw.table.Point2DKey.addFields( 

91 self.schema, "b", "", "pixel") 

92 self.xErrKey = self.schema.addField("b_xErr", type="F") 

93 self.yErrKey = self.schema.addField("b_yErr", type="F") 

94 self.centroidFlagKey = self.schema.addField("b_flag", type="Flag") 

95 

96 self.shapeKey = lsst.afw.table.QuadrupoleKey.addFields( 

97 self.schema, "c", "", lsst.afw.table.CoordinateType.PIXEL) 

98 self.xxErrKey = self.schema.addField("c_xxErr", type="F") 

99 self.xyErrKey = self.schema.addField("c_xyErr", type="F") 

100 self.yyErrKey = self.schema.addField("c_yyErr", type="F") 

101 self.shapeFlagKey = self.schema.addField("c_flag", type="Flag") 

102 

103 self.psfShapeKey = lsst.afw.table.QuadrupoleKey.addFields( 

104 self.schema, "d", "", lsst.afw.table.CoordinateType.PIXEL) 

105 self.psfShapeFlagKey = self.schema.addField("d_flag", type="Flag") 

106 

107 self.coordErrKey = lsst.afw.table.CoordKey.addErrorFields(self.schema) 

108 

109 self.table = lsst.afw.table.SourceTable.make(self.schema) 

110 self.catalog = lsst.afw.table.SourceCatalog(self.table) 

111 self.record = self.catalog.addNew() 

112 self.fillRecord(self.record) 

113 self.record.setId(50) 

114 self.fillRecord(self.catalog.addNew()) 

115 self.fillRecord(self.catalog.addNew()) 

116 

117 def tearDown(self): 

118 del self.schema 

119 del self.record 

120 del self.table 

121 del self.catalog 

122 

123 def checkCanonical(self): 

124 self.assertEqual(self.record.get(self.instFluxKey), 

125 self.record.getPsfInstFlux()) 

126 self.assertEqual(self.record.get(self.fluxFlagKey), 

127 self.record.getPsfFluxFlag()) 

128 self.assertEqual(self.table.getSchema().getAliasMap().get("slot_Centroid"), "b") 

129 self.assertEqual(self.centroidKey.get(self.record), 

130 self.record.getCentroid()) 

131 self.assertFloatsAlmostEqual( 

132 math.fabs(self.record.get(self.xErrKey)), 

133 math.sqrt(self.record.getCentroidErr()[0, 0]), rtol=1e-6) 

134 self.assertFloatsAlmostEqual( 

135 math.fabs(self.record.get(self.yErrKey)), 

136 math.sqrt(self.record.getCentroidErr()[1, 1]), rtol=1e-6) 

137 self.assertEqual(self.table.getSchema().getAliasMap().get("slot_Shape"), "c") 

138 self.assertEqual(self.shapeKey.get(self.record), 

139 self.record.getShape()) 

140 self.assertFloatsAlmostEqual( 

141 math.fabs(self.record.get(self.xxErrKey)), 

142 math.sqrt(self.record.getShapeErr()[0, 0]), rtol=1e-6) 

143 self.assertFloatsAlmostEqual( 

144 math.fabs(self.record.get(self.yyErrKey)), 

145 math.sqrt(self.record.getShapeErr()[1, 1]), rtol=1e-6) 

146 self.assertFloatsAlmostEqual( 

147 math.fabs(self.record.get(self.xyErrKey)), 

148 math.sqrt(self.record.getShapeErr()[2, 2]), rtol=1e-6) 

149 self.assertEqual(self.table.getSchema().getAliasMap().get("slot_PsfShape"), "d") 

150 self.assertEqual(self.psfShapeKey.get(self.record), 

151 self.record.getPsfShape()) 

152 

153 def testPersisted(self): 

154 self.table.definePsfFlux("a") 

155 self.table.defineCentroid("b") 

156 self.table.defineShape("c") 

157 self.table.definePsfShape("d") 

158 with lsst.utils.tests.getTempFilePath(".fits") as filename: 

159 self.catalog.writeFits(filename) 

160 catalog = lsst.afw.table.SourceCatalog.readFits(filename) 

161 table = catalog.getTable() 

162 record = catalog[0] 

163 # I'm using the keys from the non-persisted table. They should work at least in the 

164 # current implementation 

165 self.assertEqual(record.get(self.instFluxKey), record.getPsfInstFlux()) 

166 self.assertEqual(record.get(self.fluxFlagKey), record.getPsfFluxFlag()) 

167 self.assertEqual(table.getSchema().getAliasMap().get("slot_Centroid"), "b") 

168 centroid = self.centroidKey.get(self.record) 

169 self.assertEqual(centroid, record.getCentroid()) 

170 self.assertFloatsAlmostEqual( 

171 math.fabs(self.record.get(self.xErrKey)), 

172 math.sqrt(self.record.getCentroidErr()[0, 0]), rtol=1e-6) 

173 self.assertFloatsAlmostEqual( 

174 math.fabs(self.record.get(self.yErrKey)), 

175 math.sqrt(self.record.getCentroidErr()[1, 1]), rtol=1e-6) 

176 shape = self.shapeKey.get(self.record) 

177 self.assertEqual(table.getSchema().getAliasMap().get("slot_Shape"), "c") 

178 self.assertEqual(shape, record.getShape()) 

179 self.assertFloatsAlmostEqual( 

180 math.fabs(self.record.get(self.xxErrKey)), 

181 math.sqrt(self.record.getShapeErr()[0, 0]), rtol=1e-6) 

182 self.assertFloatsAlmostEqual( 

183 math.fabs(self.record.get(self.yyErrKey)), 

184 math.sqrt(self.record.getShapeErr()[1, 1]), rtol=1e-6) 

185 self.assertFloatsAlmostEqual( 

186 math.fabs(self.record.get(self.xyErrKey)), 

187 math.sqrt(self.record.getShapeErr()[2, 2]), rtol=1e-6) 

188 psfShape = self.psfShapeKey.get(self.record) 

189 self.assertEqual(table.getSchema().getAliasMap().get("slot_PsfShape"), "d") 

190 self.assertEqual(psfShape, record.getPsfShape()) 

191 

192 def testCanonical2(self): 

193 self.table.definePsfFlux("a") 

194 self.table.defineCentroid("b") 

195 self.table.defineShape("c") 

196 self.table.definePsfShape("d") 

197 self.checkCanonical() 

198 

199 def testPickle(self): 

200 p = pickle.dumps(self.catalog) 

201 new = pickle.loads(p) 

202 

203 self.assertEqual(self.catalog.schema.getNames(), new.schema.getNames()) 

204 self.assertEqual(len(self.catalog), len(new)) 

205 for r1, r2 in zip(self.catalog, new): 

206 # Columns that are easy to test 

207 for field in ("a_instFlux", "a_instFluxErr", "id"): 

208 k1 = self.catalog.schema.find(field).getKey() 

209 k2 = new.schema.find(field).getKey() 

210 self.assertEqual(r1[k1], r2[k2]) 

211 

212 def testCoordUpdate(self): 

213 self.table.defineCentroid("b") 

214 wcs = makeWcs() 

215 self.record.updateCoord(wcs) 

216 coord1 = self.record.getCoord() 

217 coord2 = wcs.pixelToSky(self.record.get(self.centroidKey)) 

218 self.assertEqual(coord1, coord2) 

219 

220 def testCoordErrors(self): 

221 self.table.defineCentroid("b") 

222 wcs = makeWcs() 

223 self.record.updateCoord(wcs) 

224 

225 scale = (1.0 * lsst.geom.arcseconds).asDegrees() 

226 center = self.record.getCentroid() 

227 skyCenter = wcs.pixelToSky(center) 

228 localGnomonicWcs = lsst.afw.geom.makeSkyWcs( 

229 center, skyCenter, np.diag((scale, scale))) 

230 measurementToLocalGnomonic = wcs.getTransform().then( 

231 localGnomonicWcs.getTransform().inverted() 

232 ) 

233 localMatrix = measurementToLocalGnomonic.getJacobian(center) 

234 radMatrix = np.radians(localMatrix / 3600) 

235 

236 centroidErr = self.record.getCentroidErr() 

237 coordErr = radMatrix.dot(centroidErr.dot(radMatrix.T)) 

238 # multiply RA by cos(Dec): 

239 cosDec = np.cos(skyCenter.getDec().asRadians()) 

240 decCorr = np.array([[cosDec, 0], [0, 1]]) 

241 coordErrCorr = decCorr.dot(coordErr.dot(decCorr)) 

242 catCoordErr = self.record.get(self.coordErrKey) 

243 np.testing.assert_almost_equal(coordErrCorr, catCoordErr, decimal=16) 

244 

245 def testSorting(self): 

246 self.assertFalse(self.catalog.isSorted()) 

247 self.catalog.sort() 

248 self.assertTrue(self.catalog.isSorted()) 

249 r = self.catalog.find(2) 

250 self.assertEqual(r["id"], 2) 

251 r = self.catalog.find(500) 

252 self.assertIsNone(r) 

253 

254 def testConversion(self): 

255 catalog1 = self.catalog.cast(lsst.afw.table.SourceCatalog) 

256 catalog2 = self.catalog.cast(lsst.afw.table.SimpleCatalog) 

257 catalog3 = self.catalog.cast(lsst.afw.table.SourceCatalog, deep=True) 

258 catalog4 = self.catalog.cast(lsst.afw.table.SimpleCatalog, deep=True) 

259 self.assertEqual(self.catalog.table, catalog1.table) 

260 self.assertEqual(self.catalog.table, catalog2.table) 

261 self.assertNotEqual(self.catalog.table, catalog3.table) 

262 self.assertNotEqual(self.catalog.table, catalog3.table) 

263 for r, r1, r2, r3, r4 in zip(self.catalog, catalog1, catalog2, catalog3, catalog4): 

264 self.assertEqual(r, r1) 

265 self.assertEqual(r, r2) 

266 self.assertNotEqual(r, r3) 

267 self.assertNotEqual(r, r4) 

268 self.assertEqual(r.getId(), r3.getId()) 

269 self.assertEqual(r.getId(), r4.getId()) 

270 

271 def testColumnView(self): 

272 cols1 = self.catalog.getColumnView() 

273 cols2 = self.catalog.columns 

274 self.assertIs(cols1, cols2) 

275 self.assertIsInstance(cols1, lsst.afw.table.SourceColumnView) 

276 self.table.definePsfFlux("a") 

277 self.table.defineCentroid("b") 

278 self.table.defineShape("c") 

279 self.table.definePsfShape("d") 

280 self.assertFloatsEqual(cols2["a_instFlux"], cols2.getPsfInstFlux()) 

281 self.assertFloatsEqual(cols2["a_instFluxErr"], cols2.getPsfInstFluxErr()) 

282 self.assertFloatsEqual(cols2["b_x"], cols2.getX()) 

283 self.assertFloatsEqual(cols2["b_y"], cols2.getY()) 

284 self.assertFloatsEqual(cols2["c_xx"], cols2.getIxx()) 

285 self.assertFloatsEqual(cols2["c_yy"], cols2.getIyy()) 

286 self.assertFloatsEqual(cols2["c_xy"], cols2.getIxy()) 

287 self.assertFloatsEqual(cols2["d_xx"], cols2.getPsfIxx()) 

288 self.assertFloatsEqual(cols2["d_yy"], cols2.getPsfIyy()) 

289 self.assertFloatsEqual(cols2["d_xy"], cols2.getPsfIxy()) 

290 

291 # Trying to access slots which have been removed should raise. 

292 self.catalog.table.schema.getAliasMap().erase("slot_Centroid") 

293 self.catalog.table.schema.getAliasMap().erase("slot_Shape") 

294 self.catalog.table.schema.getAliasMap().erase("slot_PsfShape") 

295 for quantity in ["X", "Y", "Ixx", "Iyy", "Ixy", "PsfIxx", "PsfIyy", "PsfIxy"]: 

296 with self.assertRaises(lsst.pex.exceptions.LogicError): 

297 getattr(self.catalog, f"get{quantity}")() 

298 

299 def testForwarding(self): 

300 """Verify that Catalog forwards unknown methods to its table and/or columns.""" 

301 self.table.definePsfFlux("a") 

302 self.table.defineCentroid("b") 

303 self.table.defineShape("c") 

304 self.assertFloatsEqual(self.catalog.columns["a_instFlux"], 

305 self.catalog["a_instFlux"]) 

306 self.assertFloatsEqual(self.catalog.columns[self.instFluxKey], 

307 self.catalog.get(self.instFluxKey)) 

308 self.assertFloatsEqual(self.catalog.columns.get(self.instFluxKey), 

309 self.catalog.getPsfInstFlux()) 

310 self.assertEqual(self.instFluxKey, self.catalog.getPsfFluxSlot().getMeasKey()) 

311 with self.assertRaises(AttributeError): 

312 self.catalog.foo() 

313 

314 def testBitsColumn(self): 

315 

316 allBits = self.catalog.getBits() 

317 someBits = self.catalog.getBits(["a_flag", "c_flag"]) 

318 self.assertEqual(allBits.getMask("a_flag"), 0x1) 

319 self.assertEqual(allBits.getMask("b_flag"), 0x2) 

320 self.assertEqual(allBits.getMask("c_flag"), 0x4) 

321 self.assertEqual(someBits.getMask(self.fluxFlagKey), 0x1) 

322 self.assertEqual(someBits.getMask(self.shapeFlagKey), 0x2) 

323 np.testing.assert_array_equal((allBits.array & 0x1 != 0), self.catalog["a_flag"]) 

324 np.testing.assert_array_equal((allBits.array & 0x2 != 0), self.catalog["b_flag"]) 

325 np.testing.assert_array_equal((allBits.array & 0x4 != 0), self.catalog["c_flag"]) 

326 np.testing.assert_array_equal((someBits.array & 0x1 != 0), self.catalog["a_flag"]) 

327 np.testing.assert_array_equal((someBits.array & 0x2 != 0), self.catalog["c_flag"]) 

328 

329 def testCast(self): 

330 baseCat = self.catalog.cast(lsst.afw.table.BaseCatalog) 

331 baseCat.cast(lsst.afw.table.SourceCatalog) 

332 

333 def testFootprints(self): 

334 '''Test round-tripping Footprints (inc. HeavyFootprints) to FITS 

335 ''' 

336 src1 = self.catalog.addNew() 

337 src2 = self.catalog.addNew() 

338 src3 = self.catalog.addNew() 

339 self.fillRecord(src1) 

340 self.fillRecord(src2) 

341 self.fillRecord(src3) 

342 src2.setParent(src1.getId()) 

343 

344 W, H = 100, 100 

345 mim = lsst.afw.image.MaskedImageF(W, H) 

346 im = mim.getImage() 

347 msk = mim.getMask() 

348 var = mim.getVariance() 

349 x, y = np.meshgrid(np.arange(W, dtype=int), np.arange(H, dtype=int)) 

350 im.array[:] = y*1E6 + x*1E3 

351 msk.array[:] = (y << 8) | x 

352 var.array[:] = y*1E2 + x 

353 spanSet = lsst.afw.geom.SpanSet.fromShape(20).shiftedBy(50, 50) 

354 circ = lsst.afw.detection.Footprint(spanSet) 

355 heavy = lsst.afw.detection.makeHeavyFootprint(circ, mim) 

356 src2.setFootprint(heavy) 

357 

358 for i, src in enumerate(self.catalog): 

359 if src != src2: 

360 spanSet = lsst.afw.geom.SpanSet.fromShape(1 + i*2).shiftedBy(50, 50) 

361 src.setFootprint(lsst.afw.detection.Footprint(spanSet)) 

362 

363 # insert this HeavyFootprint into an otherwise blank image (for comparing the results) 

364 mim2 = lsst.afw.image.MaskedImageF(W, H) 

365 heavy.insert(mim2) 

366 

367 with lsst.utils.tests.getTempFilePath(".fits") as fn: 

368 self.catalog.writeFits(fn) 

369 

370 cat2 = lsst.afw.table.SourceCatalog.readFits(fn) 

371 r2 = cat2[-2] 

372 h2 = r2.getFootprint() 

373 self.assertTrue(h2.isHeavy()) 

374 mim3 = lsst.afw.image.MaskedImageF(W, H) 

375 h2.insert(mim3) 

376 

377 self.assertFalse(cat2[-1].getFootprint().isHeavy()) 

378 self.assertFalse(cat2[-3].getFootprint().isHeavy()) 

379 self.assertFalse(cat2[0].getFootprint().isHeavy()) 

380 self.assertFalse(cat2[1].getFootprint().isHeavy()) 

381 self.assertFalse(cat2[2].getFootprint().isHeavy()) 

382 

383 if False: 

384 # Write out before-n-after FITS images 

385 for MI in [mim, mim2, mim3]: 

386 f, fn2 = tempfile.mkstemp(prefix='testHeavyFootprint-', suffix='.fits') 

387 os.close(f) 

388 MI.writeFits(fn2) 

389 print('wrote', fn2) 

390 

391 self.assertFloatsEqual(mim2.getImage().getArray(), mim3.getImage().getArray()) 

392 self.assertFloatsEqual(mim2.getMask().getArray(), mim3.getMask().getArray()) 

393 self.assertFloatsEqual(mim2.getVariance().getArray(), mim3.getVariance().getArray()) 

394 

395 im3 = mim3.getImage() 

396 ma3 = mim3.getMask() 

397 va3 = mim3.getVariance() 

398 for y in range(H): 

399 for x in range(W): 

400 if circ.contains(lsst.geom.Point2I(x, y)): 

401 self.assertEqual(im[x, y, lsst.afw.image.PARENT], im3[x, y, lsst.afw.image.PARENT]) 

402 self.assertEqual(msk[x, y, lsst.afw.image.PARENT], ma3[x, y, lsst.afw.image.PARENT]) 

403 self.assertEqual(var[x, y, lsst.afw.image.PARENT], va3[x, y, lsst.afw.image.PARENT]) 

404 else: 

405 self.assertEqual(im3[x, y, lsst.afw.image.PARENT], 0.) 

406 self.assertEqual(ma3[x, y, lsst.afw.image.PARENT], 0.) 

407 self.assertEqual(va3[x, y, lsst.afw.image.PARENT], 0.) 

408 

409 cat3 = lsst.afw.table.SourceCatalog.readFits( 

410 fn, flags=lsst.afw.table.SOURCE_IO_NO_HEAVY_FOOTPRINTS) 

411 for src in cat3: 

412 self.assertFalse(src.getFootprint().isHeavy()) 

413 cat4 = lsst.afw.table.SourceCatalog.readFits( 

414 fn, flags=lsst.afw.table.SOURCE_IO_NO_FOOTPRINTS) 

415 for src in cat4: 

416 self.assertIsNone(src.getFootprint()) 

417 

418 self.catalog.writeFits( 

419 fn, flags=lsst.afw.table.SOURCE_IO_NO_HEAVY_FOOTPRINTS) 

420 cat5 = lsst.afw.table.SourceCatalog.readFits(fn) 

421 for src in cat5: 

422 self.assertFalse(src.getFootprint().isHeavy()) 

423 

424 self.catalog.writeFits( 

425 fn, flags=lsst.afw.table.SOURCE_IO_NO_FOOTPRINTS) 

426 cat6 = lsst.afw.table.SourceCatalog.readFits(fn) 

427 for src in cat6: 

428 self.assertIsNone(src.getFootprint()) 

429 

430 # Insert a source with no Footprint to ensure that a None footprint 

431 # is persisted correctly 

432 self.catalog.addNew() 

433 

434 with lsst.utils.tests.getTempFilePath(".fits") as fn: 

435 self.catalog.writeFits(fn) 

436 

437 cat2 = lsst.afw.table.SourceCatalog.readFits(fn) 

438 self.assertIsNone(cat2[-1].getFootprint()) 

439 

440 def testFootprintsToNumpy(self): 

441 schema = lsst.afw.table.SourceTable.makeMinimalSchema() 

442 table = lsst.afw.table.SourceTable.make(schema) 

443 catalog = lsst.afw.table.SourceCatalog(table) 

444 src1 = catalog.addNew() 

445 src2 = catalog.addNew() 

446 src3 = catalog.addNew() 

447 

448 spans1 = lsst.afw.geom.SpanSet.fromShape(5, lsst.afw.geom.Stencil.BOX).shiftedBy(10, 20) 

449 spans2 = lsst.afw.geom.SpanSet.fromShape(5, lsst.afw.geom.Stencil.BOX).shiftedBy(30, 40) 

450 spans3 = lsst.afw.geom.SpanSet.fromShape(3, lsst.afw.geom.Stencil.BOX).shiftedBy(40, 50) 

451 

452 src1.setFootprint(lsst.afw.detection.Footprint(spans1)) 

453 src2.setFootprint(lsst.afw.detection.Footprint(spans2)) 

454 src3.setFootprint(lsst.afw.detection.Footprint(spans3)) 

455 

456 truth = np.zeros((100, 100), dtype=int) 

457 truth[15:26, 5:16] = src1.getId() 

458 truth[35:46, 25:36] = src2.getId() 

459 truth[47:54, 37:44] = src3.getId() 

460 

461 footprintIds = lsst.afw.detection.footprintsToNumpy(catalog, shape=(100, 100), asBool=False) 

462 booleanArray = lsst.afw.detection.footprintsToNumpy(catalog, shape=(100, 100), asBool=True) 

463 

464 np.testing.assert_array_equal(footprintIds, truth) 

465 np.testing.assert_array_equal(booleanArray, truth > 0) 

466 

467 def testIdFactory(self): 

468 expId = int(1257198) 

469 reserved = 32 

470 factory = lsst.afw.table.IdFactory.makeSource(expId, reserved) 

471 id1 = factory() 

472 id2 = factory() 

473 self.assertEqual(id2 - id1, 1) 

474 factory.notify(0xFFFFFFFF) 

475 with self.assertRaises(lsst.pex.exceptions.LengthError): 

476 factory() 

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

478 factory.notify(0x1FFFFFFFF) 

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

480 lsst.afw.table.IdFactory.makeSource(0x1FFFFFFFF, reserved) 

481 

482 def testFamilies(self): 

483 self.catalog.sort() 

484 parents = self.catalog.getChildren(0) 

485 self.assertEqual(list(parents), list(self.catalog)) 

486 parentKey = lsst.afw.table.SourceTable.getParentKey() 

487 for parent in parents: 

488 self.assertEqual(parent.get(parentKey), 0) 

489 for i in range(10): 

490 child = self.catalog.addNew() 

491 self.fillRecord(child) 

492 child.set(parentKey, parent.getId()) 

493 childrenIter = self.catalog.getChildren([parent.getId() for parent in parents], 

494 [record.getId() for record in self.catalog]) 

495 for parent, (children, ids) in zip(parents, childrenIter): 

496 self.assertEqual(len(children), 10) 

497 self.assertEqual(len(children), len(ids)) 

498 for child, id in zip(children, ids): 

499 self.assertEqual(child.getParent(), parent.getId()) 

500 self.assertEqual(child.getId(), id) 

501 

502 # Check detection of unsorted catalog 

503 self.catalog.sort(self.instFluxKey) 

504 with self.assertRaises(AssertionError): 

505 self.catalog.getChildren(0) 

506 self.catalog.sort(parentKey) 

507 self.catalog.getChildren(0) # Just care this succeeds 

508 

509 def testFitsReadVersion0Compatibility(self): 

510 cat = lsst.afw.table.SourceCatalog.readFits(os.path.join(testPath, "data/empty-v0.fits")) 

511 self.assertTrue(cat.getPsfFluxSlot().isValid()) 

512 self.assertTrue(cat.getApFluxSlot().isValid()) 

513 self.assertTrue(cat.getGaussianFluxSlot().isValid()) 

514 self.assertTrue(cat.getModelFluxSlot().isValid()) 

515 self.assertTrue(cat.getCentroidSlot().isValid()) 

516 self.assertTrue(cat.getShapeSlot().isValid()) 

517 self.assertEqual(cat.getPsfFluxSlot().getMeasKey(), 

518 cat.schema.find("flux_psf").key) 

519 self.assertEqual(cat.getApFluxSlot().getMeasKey(), 

520 cat.schema.find("flux_sinc").key) 

521 self.assertEqual(cat.getGaussianFluxSlot().getMeasKey(), 

522 cat.schema.find("flux_naive").key) 

523 self.assertEqual(cat.getModelFluxSlot().getMeasKey(), 

524 cat.schema.find("cmodel_flux").key) 

525 self.assertEqual(cat.getCentroidSlot().getMeasKey().getX(), 

526 cat.schema.find("centroid_sdss_x").key) 

527 self.assertEqual(cat.getCentroidSlot().getMeasKey().getY(), 

528 cat.schema.find("centroid_sdss_y").key) 

529 self.assertEqual(cat.getShapeSlot().getMeasKey().getIxx(), 

530 cat.schema.find("shape_hsm_moments_xx").key) 

531 self.assertEqual(cat.getShapeSlot().getMeasKey().getIyy(), 

532 cat.schema.find("shape_hsm_moments_yy").key) 

533 self.assertEqual(cat.getShapeSlot().getMeasKey().getIxy(), 

534 cat.schema.find("shape_hsm_moments_xy").key) 

535 self.assertEqual(cat.getPsfFluxSlot().getErrKey(), 

536 cat.schema.find("flux_psf_err").key) 

537 self.assertEqual(cat.getApFluxSlot().getErrKey(), 

538 cat.schema.find("flux_sinc_err").key) 

539 self.assertEqual(cat.getGaussianFluxSlot().getErrKey(), 

540 cat.schema.find("flux_naive_err").key) 

541 self.assertEqual(cat.getModelFluxSlot().getErrKey(), 

542 cat.schema.find("cmodel_flux_err").key) 

543 self.assertEqual( 

544 cat.getCentroidSlot().getErrKey(), 

545 lsst.afw.table.CovarianceMatrix2fKey(cat.schema["centroid_sdss_err"], ["x", "y"])) 

546 self.assertEqual( 

547 cat.getShapeSlot().getErrKey(), 

548 lsst.afw.table.CovarianceMatrix3fKey(cat.schema["shape_hsm_moments_err"], ["xx", "yy", "xy"])) 

549 self.assertEqual(cat.getPsfFluxSlot().getFlagKey(), 

550 cat.schema.find("flux_psf_flags").key) 

551 self.assertEqual(cat.getApFluxSlot().getFlagKey(), 

552 cat.schema.find("flux_sinc_flags").key) 

553 self.assertEqual(cat.getGaussianFluxSlot().getFlagKey(), 

554 cat.schema.find("flux_naive_flags").key) 

555 self.assertEqual(cat.getModelFluxSlot().getFlagKey(), 

556 cat.schema.find("cmodel_flux_flags").key) 

557 self.assertEqual(cat.getCentroidSlot().getFlagKey(), 

558 cat.schema.find("centroid_sdss_flags").key) 

559 self.assertEqual(cat.getShapeSlot().getFlagKey(), 

560 cat.schema.find("shape_hsm_moments_flags").key) 

561 

562 def testFitsReadVersion1Compatibility(self): 

563 """Test reading of catalogs with version 1 schema 

564 

565 Version 1 catalogs need to have added aliases from Sigma->Err and 

566 from `_flux`->`_instFlux`. 

567 """ 

568 cat = lsst.afw.table.SourceCatalog.readFits( 

569 os.path.join(testPath, "data", "sourceTable-v1.fits")) 

570 self.assertEqual( 

571 cat.getCentroidSlot().getErrKey(), 

572 lsst.afw.table.CovarianceMatrix2fKey( 

573 cat.schema["slot_Centroid"], 

574 ["x", "y"])) 

575 self.assertEqual( 

576 cat.getShapeSlot().getErrKey(), 

577 lsst.afw.table.CovarianceMatrix3fKey(cat.schema["slot_Shape"], ["xx", "yy", "xy"])) 

578 # check the flux->instFlux conversion 

579 self.assertEqual(cat.schema["a_flux"].asKey(), cat.schema["a_instFlux"].asKey()) 

580 self.assertEqual(cat.schema["a_fluxSigma"].asKey(), cat.schema["a_instFluxErr"].asKey()) 

581 

582 def testFitsReadVersion2CompatibilityRealSourceCatalog(self): 

583 """DM-15891: some fields were getting aliases they shouldn't have.""" 

584 cat = lsst.afw.table.SourceCatalog.readFits( 

585 os.path.join(testPath, "data", "sourceCatalog-hsc-v2.fits")) 

586 self.assertNotIn('base_SdssShape_flux_xxinstFlux', cat.schema) 

587 self.assertIn('base_SdssShape_instFlux_xx_Cov', cat.schema) 

588 self.assertNotIn('base_Blendedness_abs_flux_cinstFlux', cat.schema) 

589 

590 def testFitsReadVersion2CompatibilityRealCoaddMeasCatalog(self): 

591 """DM-16068: some fields were not getting aliases they should have 

592 

593 In particular, the alias setting relies on flux fields having their 

594 units set properly. Prior to a resolution of DM-16068, the units 

595 were not getting set for several CModel flux fields and one deblender 

596 field (deblend_psfFlux), and thus were not getting the 

597 `_flux`->`_instFlux` aliases set. 

598 

599 NOTE: this test will (is meant to) fail on the read until DM-16068 is 

600 resolved. The error is: 

601 

602 lsst::pex::exceptions::NotFoundError: 'Field or subfield with 

603 name 'modelfit_CModel_instFlux' not found with type 'D'.' 

604 """ 

605 cat = lsst.afw.table.SourceCatalog.readFits( 

606 os.path.join(testPath, "data", "deepCoadd_meas_HSC_v2.fits")) 

607 self.assertIn('modelfit_CModel_instFlux', cat.schema) 

608 self.assertIn('modelfit_CModel_instFluxErr', cat.schema) 

609 self.assertIn('modelfit_CModel_instFlux_inner', cat.schema) 

610 self.assertIn('modelfit_CModel_dev_instFlux_inner', cat.schema) 

611 self.assertIn('modelfit_CModel_exp_instFlux_inner', cat.schema) 

612 self.assertIn('modelfit_CModel_initial_instFlux_inner', cat.schema) 

613 

614 def testFitsVersion2Compatibility(self): 

615 """Test reading of catalogs with version 2 schema 

616 

617 Version 2 catalogs need to have added aliases from `_flux`->`_instFlux`. 

618 """ 

619 cat = lsst.afw.table.SourceCatalog.readFits(os.path.join(testPath, "data", "sourceTable-v2.fits")) 

620 # check the flux->instFlux conversion 

621 self.assertEqual(cat.schema["a_flux"].asKey(), cat.schema["a_instFlux"].asKey()) 

622 self.assertEqual(cat.schema["a_fluxErr"].asKey(), cat.schema["a_instFluxErr"].asKey()) 

623 

624 def testDM1083(self): 

625 schema = lsst.afw.table.SourceTable.makeMinimalSchema() 

626 st = lsst.afw.table.SourceTable.make(schema) 

627 cat = lsst.afw.table.SourceCatalog(st) 

628 tmp = lsst.afw.table.SourceCatalog(cat.getTable()) 

629 record = tmp.addNew() 

630 cat.extend(tmp) 

631 self.assertEqual(cat[0].getId(), record.getId()) 

632 # check that the same record is in both catalogs (not a copy) 

633 record.setId(15) 

634 self.assertEqual(cat[0].getId(), record.getId()) 

635 

636 def testSlotUndefine(self): 

637 """Test that we can correctly define and undefine a slot after a SourceTable has been created""" 

638 schema = lsst.afw.table.SourceTable.makeMinimalSchema() 

639 key = schema.addField("a_instFlux", type=np.float64, doc="flux field") 

640 table = lsst.afw.table.SourceTable.make(schema) 

641 table.definePsfFlux("a") 

642 self.assertEqual(table.getPsfFluxSlot().getMeasKey(), key) 

643 table.schema.getAliasMap().erase("slot_PsfFlux") 

644 self.assertFalse(table.getPsfFluxSlot().isValid()) 

645 

646 def testOldFootprintPersistence(self): 

647 """Test that we can still read SourceCatalogs with (Heavy)Footprints saved by an older 

648 version of the pipeline with a different format. 

649 """ 

650 filename = os.path.join(testPath, "data", "old-footprint-persistence.fits") 

651 catalog1 = lsst.afw.table.SourceCatalog.readFits(filename) 

652 self.assertEqual(len(catalog1), 2) 

653 with self.assertRaises(LookupError): 

654 catalog1.schema.find("footprint") 

655 fp1 = catalog1[0].getFootprint() 

656 fp2 = catalog1[1].getFootprint() 

657 self.assertEqual(fp1.getArea(), 495) 

658 self.assertEqual(fp2.getArea(), 767) 

659 self.assertFalse(fp1.isHeavy()) 

660 self.assertTrue(fp2.isHeavy()) 

661 self.assertEqual(len(fp1.getSpans()), 29) 

662 self.assertEqual(len(fp2.getSpans()), 44) 

663 self.assertEqual(len(fp1.getPeaks()), 1) 

664 self.assertEqual(len(fp2.getPeaks()), 1) 

665 self.assertEqual(fp1.getBBox(), 

666 lsst.geom.Box2I(lsst.geom.Point2I(129, 2), lsst.geom.Extent2I(25, 29))) 

667 self.assertEqual(fp2.getBBox(), 

668 lsst.geom.Box2I(lsst.geom.Point2I(1184, 2), lsst.geom.Extent2I(78, 38))) 

669 hfp = lsst.afw.detection.HeavyFootprintF(fp2) 

670 self.assertEqual(len(hfp.getImageArray()), fp2.getArea()) 

671 self.assertEqual(len(hfp.getMaskArray()), fp2.getArea()) 

672 self.assertEqual(len(hfp.getVarianceArray()), fp2.getArea()) 

673 catalog2 = lsst.afw.table.SourceCatalog.readFits( 

674 filename, flags=lsst.afw.table.SOURCE_IO_NO_HEAVY_FOOTPRINTS) 

675 self.assertEqual(list(fp1.getSpans()), 

676 list(catalog2[0].getFootprint().getSpans())) 

677 self.assertEqual(list(fp2.getSpans()), 

678 list(catalog2[1].getFootprint().getSpans())) 

679 self.assertFalse(catalog2[1].getFootprint().isHeavy()) 

680 catalog3 = lsst.afw.table.SourceCatalog.readFits( 

681 filename, flags=lsst.afw.table.SOURCE_IO_NO_FOOTPRINTS) 

682 self.assertEqual(catalog3[0].getFootprint(), None) 

683 self.assertEqual(catalog3[1].getFootprint(), None) 

684 

685 def _testFluxSlot(self, slotName): 

686 """Demonstrate that we can create & use the named Flux slot.""" 

687 schema = lsst.afw.table.SourceTable.makeMinimalSchema() 

688 baseName = "afw_Test" 

689 instFluxKey = schema.addField(f"{baseName}_instFlux", type=np.float64, doc="flux") 

690 errKey = schema.addField(f"{baseName}_instFluxErr", type=np.float64, doc="flux uncertainty") 

691 flagKey = schema.addField(f"{baseName}_flag", type="Flag", doc="flux flag") 

692 catalog = lsst.afw.table.SourceCatalog(schema) 

693 table = catalog.table 

694 

695 # Initially, the slot is undefined. 

696 self.assertFalse(getattr(table, f"get{slotName}Slot")().isValid()) 

697 

698 # After definition, it maps to the keys defined above. 

699 getattr(table, f"define{slotName}")(baseName) 

700 self.assertTrue(getattr(table, f"get{slotName}Slot")().isValid()) 

701 self.assertEqual(getattr(table, f"get{slotName}Slot")().getMeasKey(), instFluxKey) 

702 self.assertEqual(getattr(table, f"get{slotName}Slot")().getErrKey(), errKey) 

703 self.assertEqual(getattr(table, f"get{slotName}Slot")().getFlagKey(), flagKey) 

704 

705 # We should be able to retrieve arbitrary values set in records. 

706 record = catalog.addNew() 

707 instFlux, err, flag = 10.0, 1.0, False 

708 record.set(instFluxKey, instFlux) 

709 record.set(errKey, err) 

710 record.set(flagKey, flag) 

711 instFluxName = slotName.replace("Flux", "InstFlux") 

712 self.assertEqual(getattr(record, f"get{instFluxName}")(), instFlux) 

713 self.assertEqual(getattr(record, f"get{instFluxName}Err")(), err) 

714 self.assertEqual(getattr(record, f"get{slotName}Flag")(), flag) 

715 

716 # Should also be able to retrieve them as columns from the catalog 

717 self.assertEqual(getattr(catalog, f"get{instFluxName}")()[0], instFlux) 

718 self.assertEqual(getattr(catalog, f"get{instFluxName}Err")()[0], err) 

719 

720 # And we should be able to delete the slot, breaking the mapping. 

721 table.schema.getAliasMap().erase(f"slot_{slotName}") 

722 self.assertFalse(getattr(table, f"get{slotName}Slot")().isValid()) 

723 self.assertNotEqual(getattr(table, f"get{slotName}Slot")().getMeasKey(), instFluxKey) 

724 self.assertNotEqual(getattr(table, f"get{slotName}Slot")().getErrKey(), errKey) 

725 self.assertNotEqual(getattr(table, f"get{slotName}Slot")().getFlagKey(), flagKey) 

726 

727 # When the slot has been deleted, attempting to access it should 

728 # throw a LogicError. 

729 with self.assertRaises(lsst.pex.exceptions.LogicError): 

730 getattr(catalog, f"get{instFluxName}")() 

731 with self.assertRaises(lsst.pex.exceptions.LogicError): 

732 getattr(catalog, f"get{instFluxName}Err")() 

733 

734 def testFluxSlots(self): 

735 """Check that all the expected flux slots are present & correct.""" 

736 for slotName in ["ApFlux", "CalibFlux", "GaussianFlux", "ModelFlux", 

737 "PsfFlux"]: 

738 self._testFluxSlot(slotName) 

739 

740 # But, of course, we should not accept a slot which hasn't be defined. 

741 with self.assertRaises(AttributeError): 

742 self._testFluxSlot("NotExtantFlux") 

743 

744 def testStr(self): 

745 """Check that the str() produced on a catalog contains expected things.""" 

746 string = str(self.catalog) 

747 for field in ('id', 'coord_ra', 'coord_dec'): 

748 self.assertIn(field, string) 

749 

750 def testRepr(self): 

751 """Check that the repr() produced on a catalog contains expected things.""" 

752 string = repr(self.catalog) 

753 self.assertIn(str(type(self.catalog)), string) 

754 for field in ('id', 'coord_ra', 'coord_dec'): 

755 self.assertIn(field, string) 

756 

757 def testStrNonContiguous(self): 

758 """Check that str() doesn't fail on non-contiguous tables.""" 

759 del self.catalog[1] 

760 string = str(self.catalog) 

761 self.assertIn('Non-contiguous afw.Catalog of 2 rows.', string) 

762 for field in ('id', 'coord_ra', 'coord_dec'): 

763 self.assertIn(field, string) 

764 

765 def testRecordStr(self): 

766 """Test that str(record) contains expected things.""" 

767 string = str(self.catalog[0]) 

768 for field in ('id: 50', 'coord_ra: nan', 'coord_dec: nan'): 

769 self.assertIn(field, string) 

770 

771 def testRecordRepr(self): 

772 """Test that repr(record) contains expected things.""" 

773 string = repr(self.catalog[0]) 

774 self.assertIn(str(type(self.catalog[0])), string) 

775 for field in ('id: 50', 'coord_ra: nan', 'coord_dec: nan'): 

776 self.assertIn(field, string) 

777 

778 def testGetNonContiguous(self): 

779 """Check that we can index on non-contiguous tables""" 

780 # Make a non-contiguous catalog 

781 nonContiguous = type(self.catalog)(self.catalog.table) 

782 for rr in reversed(self.catalog): 

783 nonContiguous.append(rr) 

784 num = len(self.catalog) 

785 # Check assumptions 

786 self.assertFalse(nonContiguous.isContiguous()) # We managed to produce a non-contiguous catalog 

787 self.assertEqual(len(set(self.catalog["id"])), num) # ID values are unique 

788 # Indexing with boolean array 

789 select = np.zeros(num, dtype=bool) 

790 select[1] = True 

791 self.assertEqual(nonContiguous[np.flip(select, 0)]["id"], self.catalog[select]["id"]) 

792 # Extracting a number column 

793 column = "a_instFlux" 

794 array = nonContiguous[column] 

795 self.assertFloatsEqual(np.flip(array, 0), self.catalog[column]) 

796 with self.assertRaises(ValueError): 

797 array[1] = 1.2345 # Should be immutable 

798 # Extracting a flag column 

799 column = "a_flag" 

800 array = nonContiguous[column] 

801 np.testing.assert_equal(np.flip(array, 0), self.catalog[column]) 

802 with self.assertRaises(ValueError): 

803 array[1] = True # Should be immutable 

804 

805 

806class MemoryTester(lsst.utils.tests.MemoryTestCase): 

807 pass 

808 

809 

810def setup_module(module): 

811 lsst.utils.tests.init() 

812 

813 

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

815 lsst.utils.tests.init() 

816 unittest.main()