Coverage for tests/test_ndf_output_archive.py: 99%

1# This file is part of lsst-images. 

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# Use of this source code is governed by a 3-clause BSD-style 

10# license that can be found in the LICENSE file. 

11 

12from __future__ import annotations 

13 

14import json 

15import tempfile 

16import unittest 

17from unittest import mock 

18 

19import astropy.io.fits 

20import astropy.table 

21import astropy.units as u 

22import numpy as np 

23import pydantic 

24 

25from lsst.images import Box, Image, MaskedImage, MaskPlane, MaskSchema 

26from lsst.images._transforms import FrameLookupError, FrameSet, Transform 

27from lsst.images._transforms._frames import DetectorFrame, Frame 

28from lsst.images.fits import ExtensionKey, FitsOpaqueMetadata 

29from lsst.images.serialization import ( 

30 ArrayReferenceModel, 

31 InlineArrayModel, 

32 open_archive, 

33 read_archive, 

34) 

35from lsst.images.tests import make_random_sky_projection 

36 

37try: 

38 import h5py 

39 

40 from lsst.images.ndf import ( 

41 NdfInputArchive, 

42 NdfOutputArchive, 

43 _hds, 

44 write, 

45 ) 

46 from lsst.images.ndf._hds import DAT__SZNAM 

47 

48 HAVE_H5PY = True 

49except ImportError: 

50 HAVE_H5PY = False 

51 

52 

53class TinyFrameSet(FrameSet): 

54 """Minimal concrete frame-set for archive bookkeeping tests.""" 

55 

56 def __contains__(self, frame: Frame) -> bool: 

57 return False 

58 

59 def __getitem__[I: Frame, O: Frame](self, key: tuple[I, O]) -> Transform[I, O]: 

60 raise FrameLookupError(key) 

61 

62 

63class TinyTree(pydantic.BaseModel): 

64 """A trivial Pydantic model used as a serialization stand-in.""" 

65 

66 name: str 

67 

68 

69@unittest.skipUnless(HAVE_H5PY, "h5py is not installed") 

70class NdfOutputArchiveBasicsTestCase(unittest.TestCase): 

71 """Tests for `NdfOutputArchive` constructor and `serialize_direct`.""" 

72 

73 def test_serialize_direct_calls_serializer_with_nested_archive(self): 

74 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

75 with h5py.File(tmp.name, "w") as f: 

76 arch = NdfOutputArchive(f) 

77 tree = arch.serialize_direct("top", lambda nested: TinyTree(name="hello")) 

78 self.assertEqual(tree.name, "hello") 

79 

80 def test_constructor_marks_root_as_ndf(self): 

81 """The constructor should set CLASS=NDF on the root group so that 

82 Starlink tools recognise the file as an NDF. 

83 """ 

84 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

85 with h5py.File(tmp.name, "w") as f: 

86 NdfOutputArchive(f) 

87 with h5py.File(tmp.name, "r") as f: 

88 self.assertEqual(f["/"].attrs["CLASS"], b"NDF") 

89 

90 

91@unittest.skipUnless(HAVE_H5PY, "h5py is not installed") 

92class NdfOutputArchiveAddArrayTestCase(unittest.TestCase): 

93 """Tests for `NdfOutputArchive.add_array` routing.""" 

94 

95 def test_top_level_image_routes_to_data_array(self): 

96 data = np.arange(20, dtype=np.float32).reshape(4, 5) 

97 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

98 with h5py.File(tmp.name, "w") as f: 

99 arch = NdfOutputArchive(f) 

100 ref = arch.add_array(data, name="image") 

101 self.assertEqual(ref.source, "ndf:/DATA_ARRAY/DATA") 

102 with h5py.File(tmp.name, "r") as f: 

103 ds = f["/DATA_ARRAY/DATA"] 

104 self.assertEqual(ds.dtype, np.float32) 

105 np.testing.assert_array_equal(ds[()], data) 

106 self.assertEqual(f["/DATA_ARRAY"].attrs["CLASS"], b"ARRAY") 

107 origin = f["/DATA_ARRAY/ORIGIN"] 

108 self.assertEqual(origin.dtype, np.int64) 

109 self.assertEqual(origin.shape, (2,)) 

110 

111 def test_top_level_variance_routes_to_variance(self): 

112 data = np.full((3, 3), 0.5, dtype=np.float64) 

113 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

114 with h5py.File(tmp.name, "w") as f: 

115 arch = NdfOutputArchive(f) 

116 ref = arch.add_array(data, name="variance") 

117 self.assertEqual(ref.source, "ndf:/VARIANCE/DATA") 

118 with h5py.File(tmp.name, "r") as f: 

119 self.assertEqual(f["/VARIANCE"].attrs["CLASS"], b"ARRAY") 

120 self.assertEqual(f["/VARIANCE/DATA"].dtype, np.float64) 

121 

122 def test_top_level_compatible_mask_routes_to_quality(self): 

123 data = np.array([[0, 1, 2], [3, 4, 5]], dtype=np.uint8) 

124 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

125 with h5py.File(tmp.name, "w") as f: 

126 arch = NdfOutputArchive(f) 

127 ref = arch.add_array(data, name="mask") 

128 self.assertEqual(ref.source, "ndf:/QUALITY/QUALITY/DATA") 

129 with h5py.File(tmp.name, "r") as f: 

130 self.assertEqual(f["/QUALITY"].attrs["CLASS"], b"QUALITY") 

131 self.assertEqual(f["/QUALITY/QUALITY"].attrs["CLASS"], b"ARRAY") 

132 self.assertEqual(f["/QUALITY/QUALITY/DATA"].dtype, np.uint8) 

133 np.testing.assert_array_equal(f["/QUALITY/QUALITY/DATA"][()], data) 

134 self.assertEqual(f["/QUALITY/QUALITY/ORIGIN"].dtype, np.int32) 

135 self.assertEqual(f["/QUALITY/QUALITY/ORIGIN"].shape, (2,)) 

136 self.assertEqual(f["/QUALITY/QUALITY/BAD_PIXEL"].id.get_type().get_class(), h5py.h5t.BITFIELD) 

137 self.assertFalse(_hds.read_array(f["/QUALITY/QUALITY/BAD_PIXEL"])) 

138 self.assertEqual(f["/QUALITY/BADBITS"][()], 255) 

139 

140 def test_top_level_incompatible_mask_routes_to_more_lsst(self): 

141 # 3D mask array in NDF storage order (mask-byte, y, x) is hoisted 

142 # as a sub-NDF inside /MORE/LSST/MASK, with a compressed 2D view 

143 # exposed as /QUALITY/QUALITY for standard NDF applications. 

144 data = np.zeros((2, 3, 4), dtype=np.uint8) 

145 data[0, 1, 2] = 4 

146 data[1, 2, 3] = 8 

147 expected_quality = np.any(data != 0, axis=0).astype(np.uint8) 

148 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

149 with h5py.File(tmp.name, "w") as f: 

150 arch = NdfOutputArchive(f) 

151 ref = arch.add_array(data, name="mask") 

152 self.assertEqual(ref.source, "ndf:/MORE/LSST/MASK/DATA_ARRAY/DATA") 

153 with h5py.File(tmp.name, "r") as f: 

154 # /MORE/LSST/MASK is a real NDF: top-level CLASS="NDF" 

155 # containing a DATA_ARRAY structure with DATA + ORIGIN. 

156 self.assertEqual(f["/MORE/LSST/MASK"].attrs["CLASS"], b"NDF") 

157 self.assertEqual(f["/MORE/LSST/MASK/DATA_ARRAY"].attrs["CLASS"], b"ARRAY") 

158 self.assertEqual(f["/MORE/LSST/MASK/DATA_ARRAY/DATA"].shape, data.shape) 

159 self.assertEqual(f["/QUALITY/QUALITY"].attrs["CLASS"], b"ARRAY") 

160 np.testing.assert_array_equal(f["/QUALITY/QUALITY/DATA"][()], expected_quality) 

161 self.assertEqual(f["/QUALITY/BADBITS"][()], 255) 

162 origin = f["/MORE/LSST/MASK/DATA_ARRAY/ORIGIN"] 

163 self.assertEqual(origin.dtype, np.int64) 

164 self.assertEqual(origin.shape, (3,)) 

165 

166 def test_long_hoisted_component_is_shrunk(self): 

167 # Regression for the cell_coadd failure: the /noise_realizations/0 

168 # archive path contains an 18-character component. 

169 data = np.array([[1.0, 2.0]], dtype=np.float32) 

170 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

171 with h5py.File(tmp.name, "w") as f: 

172 arch = NdfOutputArchive(f) 

173 ref = arch.add_array(data, name="noise_realizations/0") 

174 # The reported path is exactly what is stored in the JSON. 

175 self.assertTrue(ref.source.startswith("ndf:/MORE/LSST/")) 

176 self.assertTrue(ref.source.endswith("/DATA_ARRAY/DATA")) 

177 with h5py.File(tmp.name, "r") as f: 

178 # Every HDS component is within the limit. 

179 hdf5_path = ref.source[len("ndf:") :] 

180 for component in hdf5_path.strip("/").split("/"): 

181 self.assertLessEqual(len(component), DAT__SZNAM) 

182 # The node the JSON points at actually exists. 

183 self.assertIn(hdf5_path, f) 

184 

185 def test_long_name_round_trips_through_input_archive(self): 

186 from lsst.images.ndf import NdfInputArchive 

187 

188 data = np.arange(6, dtype=np.float32).reshape(2, 3) 

189 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

190 with h5py.File(tmp.name, "w") as f: 

191 arch = NdfOutputArchive(f) 

192 ref = arch.add_array(data, name="noise_realizations/0") 

193 with NdfInputArchive.open(tmp.name) as inp: 

194 read_back = inp.get_array(ref) 

195 np.testing.assert_array_equal(read_back, data) 

196 

197 def test_repeated_long_name_gets_distinct_versioned_paths(self): 

198 data = np.array([[1.0]], dtype=np.float32) 

199 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

200 with h5py.File(tmp.name, "w") as f: 

201 arch = NdfOutputArchive(f) 

202 first = arch.add_array(data, name="noise_realizations_value") 

203 second = arch.add_array(data, name="noise_realizations_value") 

204 self.assertNotEqual(first.source, second.source) 

205 # The second occurrence keeps a visible _2 version suffix. 

206 second_leaf = second.source[len("ndf:") :].split("/")[-3] 

207 self.assertTrue(second_leaf.endswith("_2")) 

208 with h5py.File(tmp.name, "r") as f: 

209 self.assertIn(first.source[len("ndf:") :], f) 

210 self.assertIn(second.source[len("ndf:") :], f) 

211 

212 def test_nested_array_hoists_as_sub_ndf(self): 

213 # Hoisted numeric arrays land under /MORE/LSST as hierarchical 

214 # sub-NDFs (CLASS="NDF" with DATA_ARRAY/DATA + ORIGIN inside) so 

215 # Starlink tools can inspect them as ordinary NDFs while each HDS 

216 # component stays short. 

217 data = np.array([[1.0, 2.0]], dtype=np.float32) 

218 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

219 with h5py.File(tmp.name, "w") as f: 

220 arch = NdfOutputArchive(f) 

221 ref = arch.add_array(data, name="psf/coefficients") 

222 self.assertEqual(ref.source, "ndf:/MORE/LSST/PSF/COEFFICIENTS/DATA_ARRAY/DATA") 

223 with h5py.File(tmp.name, "r") as f: 

224 self.assertIn("MORE", f) 

225 self.assertIn("LSST", f["/MORE"]) 

226 self.assertIn("PSF", f["/MORE/LSST"]) 

227 self.assertIn("COEFFICIENTS", f["/MORE/LSST/PSF"]) 

228 sub = f["/MORE/LSST/PSF/COEFFICIENTS"] 

229 self.assertEqual(sub.attrs["CLASS"], b"NDF") 

230 self.assertEqual(sub["DATA_ARRAY"].attrs["CLASS"], b"ARRAY") 

231 np.testing.assert_array_equal(sub["DATA_ARRAY/DATA"][()], data) 

232 origin = sub["DATA_ARRAY/ORIGIN"] 

233 self.assertEqual(origin.dtype, np.int64) 

234 self.assertEqual(origin.shape, (data.ndim,)) 

235 

236 def test_colliding_shrunk_names_raise(self): 

237 data = np.array([[1.0]], dtype=np.float32) 

238 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

239 with h5py.File(tmp.name, "w") as f: 

240 arch = NdfOutputArchive(f) 

241 # Force both long names to shrink to the same HDS token. 

242 with mock.patch.object( 

243 arch._name_shrinker, 

244 "shrink", 

245 side_effect=lambda name, *a, **k: name.upper() if len(name) <= DAT__SZNAM else "CLASH", 

246 ): 

247 arch.add_array(data, name="long_component_name_one") 

248 with self.assertRaisesRegex(ValueError, "name collision"): 

249 arch.add_array(data, name="long_component_name_two") 

250 

251 

252@unittest.skipUnless(HAVE_H5PY, "h5py is not installed") 

253class NdfOutputArchivePointerTestCase(unittest.TestCase): 

254 """Tests for `NdfOutputArchive.serialize_pointer` and 

255 `serialize_frame_set`. 

256 """ 

257 

258 def test_serialize_pointer_writes_subtree_and_returns_pointer(self): 

259 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

260 with h5py.File(tmp.name, "w") as f: 

261 arch = NdfOutputArchive(f) 

262 ptr = arch.serialize_pointer( 

263 "psf", 

264 lambda nested: TinyTree(name="gaussian"), 

265 key=("psf", 1), 

266 ) 

267 self.assertEqual(ptr.path, "/MORE/LSST/PSF/JSON") 

268 with h5py.File(tmp.name, "r") as f: 

269 # The hoisted sub-tree is stored as a "JSON" _CHAR*N 

270 # child of the target structure. 

271 raw = f["/MORE/LSST/PSF/JSON"][()] 

272 joined = b"".join(raw).decode("ascii").rstrip(" ") 

273 self.assertIn('"name":"gaussian"', joined.replace(" ", "")) 

274 

275 def test_serialize_pointer_caches_by_key(self): 

276 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

277 with h5py.File(tmp.name, "w") as f: 

278 arch = NdfOutputArchive(f) 

279 ptr1 = arch.serialize_pointer( 

280 "psf", 

281 lambda nested: TinyTree(name="first"), 

282 key=("psf", 1), 

283 ) 

284 # Same key -> returns cached pointer; serializer not re-run 

285 # (we'd otherwise overwrite the file content with "second"). 

286 ptr2 = arch.serialize_pointer( 

287 "psf", 

288 lambda nested: TinyTree(name="second"), 

289 key=("psf", 1), 

290 ) 

291 self.assertEqual(ptr1, ptr2) 

292 with h5py.File(tmp.name, "r") as f: 

293 raw = f["/MORE/LSST/PSF/JSON"][()] 

294 joined = b"".join(raw).decode("ascii").rstrip(" ") 

295 self.assertIn("first", joined) 

296 self.assertNotIn("second", joined) 

297 

298 def test_serialize_pointer_preserves_nested_arrays(self): 

299 # Regression test: a pointer target that writes a nested array via 

300 # the nested archive must round-trip with that array still in the 

301 # file. Previously the pointer JSON was written at the target path 

302 # itself, clobbering any nested data the serializer produced. 

303 class TreeWithArray(pydantic.BaseModel): 

304 name: str 

305 data: ArrayReferenceModel 

306 

307 payload = np.arange(6, dtype=np.float32).reshape(2, 3) 

308 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

309 with h5py.File(tmp.name, "w") as f: 

310 arch = NdfOutputArchive(f) 

311 ptr = arch.serialize_pointer( 

312 "psf", 

313 lambda nested: TreeWithArray( 

314 name="gaussian", 

315 data=nested.add_array(payload, name="parameters"), 

316 ), 

317 key=("psf", 1), 

318 ) 

319 self.assertEqual(ptr.path, "/MORE/LSST/PSF/JSON") 

320 with h5py.File(tmp.name, "r") as f: 

321 # JSON stayed at <path>/JSON; the nested array is still 

322 # accessible at the sub-NDF path the serializer wrote. 

323 self.assertIn("/MORE/LSST/PSF/JSON", f) 

324 self.assertIn("/MORE/LSST/PSF/PARAMETERS/DATA_ARRAY/DATA", f) 

325 np.testing.assert_array_equal(f["/MORE/LSST/PSF/PARAMETERS/DATA_ARRAY/DATA"][()], payload) 

326 # The pointer-target structure is typed after its leaf 

327 # name, not the generic EXT. 

328 self.assertEqual(f["/MORE/LSST/PSF"].attrs["CLASS"], b"PSF") 

329 self.assertEqual(f["/MORE/LSST"].attrs["CLASS"], b"LSST") 

330 

331 def test_serialize_frame_set_records_for_iter(self): 

332 # serialize_frame_set is delegated to serialize_pointer plus 

333 # recording the (FrameSet, pointer) pair for iter_frame_sets, 

334 # mirroring how FITS and JSON archives behave. 

335 frame_set = TinyFrameSet() 

336 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

337 with h5py.File(tmp.name, "w") as f: 

338 arch = NdfOutputArchive(f) 

339 ptr = arch.serialize_frame_set( 

340 "wcs/pixel_to_sky", 

341 frame_set, 

342 lambda nested: TinyTree(name="proj"), 

343 key=("frame_set", 1), 

344 ) 

345 self.assertEqual(ptr.path, "/MORE/LSST/WCS/PIXEL_TO_SKY/JSON") 

346 recorded = list(arch.iter_frame_sets()) 

347 self.assertEqual(len(recorded), 1) 

348 self.assertIs(recorded[0][0], frame_set) 

349 self.assertEqual(recorded[0][1].path, "/MORE/LSST/WCS/PIXEL_TO_SKY/JSON") 

350 

351 

352@unittest.skipUnless(HAVE_H5PY, "h5py is not installed") 

353class NdfOutputArchiveAddTableTestCase(unittest.TestCase): 

354 """Tests for `NdfOutputArchive.add_table` and `add_structured_array`.""" 

355 

356 def test_add_table_returns_inline_table_model(self): 

357 t = astropy.table.Table({"a": [1, 2, 3], "b": [4.0, 5.0, 6.0]}) 

358 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

359 with h5py.File(tmp.name, "w") as f: 

360 arch = NdfOutputArchive(f) 

361 model = arch.add_table(t, name="some_table") 

362 self.assertEqual(len(model.columns), 2) 

363 # v1 stores tables inline in the JSON tree. 

364 self.assertIsInstance(model.columns[0].data, InlineArrayModel) 

365 

366 def test_add_structured_array_writes_column_ndfs_with_units(self): 

367 rec = np.zeros(3, dtype=[("x", np.float64), ("y", np.int32)]) 

368 rec["x"] = [1.0, 2.0, 3.0] 

369 rec["y"] = [10, 20, 30] 

370 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

371 with h5py.File(tmp.name, "w") as f: 

372 arch = NdfOutputArchive(f) 

373 model = arch.add_structured_array( 

374 rec, 

375 name="rec", 

376 units={"x": u.m}, 

377 descriptions={"y": "the y values"}, 

378 ) 

379 self.assertEqual(len(model.columns), 2) 

380 self.assertIsInstance(model.columns[0].data, ArrayReferenceModel) 

381 # Confirm units/descriptions were applied. 

382 col_x = next(c for c in model.columns if c.name == "x") 

383 col_y = next(c for c in model.columns if c.name == "y") 

384 self.assertEqual(col_x.unit, u.m) 

385 self.assertEqual(col_y.description, "the y values") 

386 self.assertEqual(col_x.data.source, "ndf:/MORE/LSST/REC/X/DATA_ARRAY/DATA") 

387 self.assertEqual(col_y.data.source, "ndf:/MORE/LSST/REC/Y/DATA_ARRAY/DATA") 

388 with h5py.File(tmp.name, "r") as f: 

389 self.assertEqual(f["/MORE/LSST/REC/X"].attrs["CLASS"], b"NDF") 

390 np.testing.assert_array_equal(f["/MORE/LSST/REC/X/DATA_ARRAY/DATA"][()], rec["x"]) 

391 self.assertEqual(f["/MORE/LSST/REC/Y"].attrs["CLASS"], b"NDF") 

392 np.testing.assert_array_equal(f["/MORE/LSST/REC/Y/DATA_ARRAY/DATA"][()], rec["y"]) 

393 with NdfInputArchive.open(tmp.name) as archive: 

394 recovered = archive.get_structured_array(model) 

395 np.testing.assert_array_equal(recovered, rec) 

396 

397 def test_add_single_column_structured_array_uses_table_name(self): 

398 rec = np.zeros(1, dtype=[("solution", np.float64, (4,))]) 

399 rec["solution"] = [[1.0, 2.0, 3.0, 4.0]] 

400 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

401 with h5py.File(tmp.name, "w") as f: 

402 arch = NdfOutputArchive(f) 

403 model = arch.add_structured_array(rec, name="psf/piff/interp/solution") 

404 self.assertEqual(len(model.columns), 1) 

405 column = model.columns[0] 

406 self.assertIsInstance(column.data, ArrayReferenceModel) 

407 self.assertEqual( 

408 column.data.source, 

409 "ndf:/MORE/LSST/PSF/PIFF/INTERP/SOLUTION/DATA_ARRAY/DATA", 

410 ) 

411 self.assertEqual(column.data.shape, [4]) 

412 with h5py.File(tmp.name, "r") as f: 

413 self.assertIn("PSF", f["/MORE/LSST"]) 

414 self.assertIn("PIFF", f["/MORE/LSST/PSF"]) 

415 self.assertIn("INTERP", f["/MORE/LSST/PSF/PIFF"]) 

416 self.assertIn("SOLUTION", f["/MORE/LSST/PSF/PIFF/INTERP"]) 

417 np.testing.assert_array_equal( 

418 f["/MORE/LSST/PSF/PIFF/INTERP/SOLUTION/DATA_ARRAY/DATA"][()], 

419 rec["solution"], 

420 ) 

421 

422 def test_structured_array_long_name_is_shrunk_and_versioned(self): 

423 dtype = np.dtype([("alpha", "f8"), ("beta", "i4")]) 

424 arr = np.zeros(3, dtype=dtype) 

425 with tempfile.NamedTemporaryFile(suffix=".sdf") as tmp: 

426 with h5py.File(tmp.name, "w") as f: 

427 arch = NdfOutputArchive(f) 

428 first = arch.add_structured_array(arr, name="catalog_of_long_named_sources") 

429 second = arch.add_structured_array(arr, name="catalog_of_long_named_sources") 

430 for model in (first, second): 

431 for column in model.columns: 

432 token = column.data.source[len("ndf:") :] 

433 for component in token.strip("/").split("/"): 

434 self.assertLessEqual(len(component), DAT__SZNAM) 

435 # The two structured arrays land in distinct sub-trees. 

436 self.assertNotEqual( 

437 first.columns[0].data.source, 

438 second.columns[0].data.source, 

439 ) 

440 # The second structured array's parent token keeps a 

441 # visible _2 version suffix. 

442 second_parent = second.columns[0].data.source[len("ndf:") :].strip("/").split("/")[-4] 

443 self.assertTrue(second_parent.endswith("_2")) 

444 with h5py.File(tmp.name, "r") as f: 

445 for model in (first, second): 

446 for column in model.columns: 

447 self.assertIn(column.data.source[len("ndf:") :], f) 

448 

449 

450@unittest.skipUnless(HAVE_H5PY, "h5py is not installed") 

451class NdfWriteWcsTestCase(unittest.TestCase): 

452 """Tests for /WCS/DATA serialization in ndf.write().""" 

453 

454 def test_write_with_projection_creates_wcs_component(self): 

455 rng = np.random.default_rng(42) 

456 det_frame = DetectorFrame(instrument="TestInst", detector=4, bbox=Box.factory[1:4096, 1:4096]) 

457 bbox = Box.factory[10:14, 20:25] 

458 sky_projection = make_random_sky_projection(rng, det_frame, Box.factory[1:4096, 1:4096]) 

459 image = Image( 

460 np.arange(20, dtype=np.float32).reshape(4, 5), 

461 bbox=bbox, 

462 sky_projection=sky_projection, 

463 ) 

464 with tempfile.NamedTemporaryFile(suffix=".sdf", delete_on_close=False) as tmp: 

465 tmp.close() 

466 write(image, tmp.name) 

467 with h5py.File(tmp.name, "r") as f: 

468 self.assertIn("WCS", f) 

469 self.assertEqual(f["/WCS"].attrs["CLASS"], b"WCS") 

470 wcs_data = f["/WCS/DATA"] 

471 self.assertEqual(wcs_data.dtype, np.dtype("|S32")) 

472 records = [s.decode("ascii").rstrip(" ") for s in wcs_data[()]] 

473 self.assertTrue(all(record[0] in {" ", "+"} for record in records)) 

474 self.assertFalse(any(record.startswith("#") for record in records)) 

475 text = _hds.decode_ndf_ast_data(records) 

476 stripped = [line.lstrip() for line in text.splitlines()] 

477 self.assertTrue(any(s.startswith("Begin FrameSet") for s in stripped)) 

478 self.assertTrue(any(s.startswith("End FrameSet") for s in stripped)) 

479 self.assertIn('Domain = "GRID"', stripped) 

480 self.assertIn('Domain = "PIXEL"', stripped) 

481 self.assertIn("Sft1 = -19", stripped) 

482 self.assertIn("Sft2 = -9", stripped) 

483 

484 def test_write_without_projection_omits_wcs_component(self): 

485 # Image with no sky_projection -> no /WCS in the file. 

486 image = Image(np.zeros((2, 2), dtype=np.float32)) 

487 with tempfile.NamedTemporaryFile(suffix=".sdf", delete_on_close=False) as tmp: 

488 tmp.close() 

489 write(image, tmp.name) 

490 with h5py.File(tmp.name, "r") as f: 

491 self.assertNotIn("WCS", f) 

492 

493 def test_mask_sub_ndf_gets_3d_wcs(self): 

494 # When an incompatible mask is hoisted to /MORE/LSST/MASK as a 

495 # sub-NDF, it should carry its own 3D /WCS. The first two axes 

496 # retain the parent image sky_projection while the third axis is 

497 # a generic mask-byte coordinate. 

498 rng = np.random.default_rng(42) 

499 det_frame = DetectorFrame(instrument="TestInst", detector=4, bbox=Box.factory[1:4096, 1:4096]) 

500 bbox = Box.factory[10:14, 20:25] 

501 sky_projection = make_random_sky_projection(rng, det_frame, Box.factory[1:4096, 1:4096]) 

502 # 12-plane schema -> native 3D uint8 mask, hoisted to /MORE/LSST/MASK. 

503 planes = [MaskPlane(f"P{i}", f"Plane {i}") for i in range(12)] 

504 image = Image( 

505 np.arange(20, dtype=np.float32).reshape(4, 5), 

506 bbox=bbox, 

507 sky_projection=sky_projection, 

508 ) 

509 masked = MaskedImage(image, mask_schema=MaskSchema(planes)) 

510 with tempfile.NamedTemporaryFile(suffix=".sdf", delete_on_close=False) as tmp: 

511 tmp.close() 

512 write(masked, tmp.name) 

513 with h5py.File(tmp.name, "r") as f: 

514 # Top-level WCS is present (existing behaviour). 

515 self.assertIn("WCS", f) 

516 top_lines = [s.decode("ascii") for s in f["/WCS/DATA"][()]] 

517 self.assertIn("MASK", f["/MORE/LSST"]) 

518 self.assertIn("WCS", f["/MORE/LSST/MASK"]) 

519 self.assertEqual(f["/MORE/LSST/MASK/WCS"].attrs["CLASS"], b"WCS") 

520 mask_lines = [s.decode("ascii") for s in f["/MORE/LSST/MASK/WCS/DATA"][()]] 

521 self.assertNotEqual(top_lines, mask_lines) 

522 mask_text = _hds.decode_ndf_ast_data(mask_lines) 

523 stripped = [line.lstrip() for line in mask_text.splitlines()] 

524 self.assertIn("Naxes = 3", stripped) 

525 self.assertIn('Domain = "GRID"', stripped) 

526 self.assertIn('Domain = "PIXEL"', stripped) 

527 self.assertIn("Sft1 = -19", stripped) 

528 self.assertIn("Sft2 = -9", stripped) 

529 self.assertIn("Sft3 = 1", stripped) 

530 self.assertIn("Begin CmpFrame", stripped) 

531 self.assertIn("Begin SkyFrame", stripped) 

532 self.assertIn('Domain = "MASK"', stripped) 

533 self.assertIn("Begin CmpMap", stripped) 

534 self.assertIn("Series = 0", stripped) 

535 

536 def test_mask_sub_ndf_no_wcs_when_image_has_no_projection(self): 

537 planes = [MaskPlane(f"P{i}", f"Plane {i}") for i in range(12)] 

538 masked = MaskedImage( 

539 Image(np.zeros((4, 5), dtype=np.float32)), 

540 mask_schema=MaskSchema(planes), 

541 ) 

542 with tempfile.NamedTemporaryFile(suffix=".sdf", delete_on_close=False) as tmp: 

543 tmp.close() 

544 write(masked, tmp.name) 

545 with h5py.File(tmp.name, "r") as f: 

546 self.assertNotIn("WCS", f) 

547 self.assertIn("MASK", f["/MORE/LSST"]) 

548 self.assertNotIn("WCS", f["/MORE/LSST/MASK"]) 

549 

550 

551@unittest.skipUnless(HAVE_H5PY, "h5py is not installed") 

552class NdfWriteFunctionTestCase(unittest.TestCase): 

553 """End-to-end tests for the module-level `write()` function.""" 

554 

555 def test_write_image_produces_valid_layout(self): 

556 image = Image( 

557 np.arange(20, dtype=np.float32).reshape(4, 5), 

558 bbox=Box.factory[10:14, 20:25], 

559 ) 

560 with tempfile.NamedTemporaryFile(suffix=".sdf", delete_on_close=False) as tmp: 

561 tmp.close() 

562 tree = write(image, tmp.name) 

563 self.assertIsNotNone(tree) 

564 with h5py.File(tmp.name, "r") as f: 

565 # Root is an NDF with a name. 

566 self.assertEqual(f["/"].attrs["CLASS"], b"NDF") 

567 self.assertIn("HDS_ROOT_NAME", f["/"].attrs) 

568 # DATA_ARRAY uses the complex form (DATA + ORIGIN). 

569 self.assertEqual(f["/DATA_ARRAY"].attrs["CLASS"], b"ARRAY") 

570 np.testing.assert_array_equal(f["/DATA_ARRAY/DATA"][()], image.array) 

571 origin = f["/DATA_ARRAY/ORIGIN"][()] 

572 self.assertEqual(origin.dtype, np.int64) 

573 self.assertEqual(len(origin), 2) 

574 # ORIGIN encodes bbox lower bounds in Fortran order. The exact 

575 # values depend on Box's API; just verify it isn't the 

576 # all-zeros placeholder when the bbox is non-trivial. 

577 self.assertFalse((origin == 0).all()) 

578 # Main JSON tree at /MORE/LSST/JSON. 

579 self.assertIn("MORE", f) 

580 self.assertIn("LSST", f["/MORE"]) 

581 self.assertIn("JSON", f["/MORE/LSST"]) 

582 

583 def test_write_image_preserves_opaque_fits_metadata(self): 

584 image = Image(np.zeros((2, 2), dtype=np.float32)) 

585 # Attach an opaque-metadata primary header to the image. 

586 primary = astropy.io.fits.Header() 

587 primary["FOO"] = ("bar", "test card") 

588 long_value = "x" * 100 

589 primary["LONGSTR"] = (long_value, "long string value") 

590 opaque = FitsOpaqueMetadata() 

591 opaque.add_header(primary, name="", ver=1) 

592 image._opaque_metadata = opaque 

593 with tempfile.NamedTemporaryFile(suffix=".sdf", delete_on_close=False) as tmp: 

594 tmp.close() 

595 write(image, tmp.name)