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 

15from pathlib import Path 

16from unittest import mock 

17 

18import astropy.io.fits 

19import astropy.table 

20import astropy.units as u 

21import numpy as np 

22import pydantic 

23import pytest 

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 

52skip_no_h5py = pytest.mark.skipif(not HAVE_H5PY, reason="h5py is not installed") 

53 

54 

55class TinyFrameSet(FrameSet): 

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

57 

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

59 return False 

60 

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

62 raise FrameLookupError(key) 

63 

64 

65class TinyTree(pydantic.BaseModel): 

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

67 

68 name: str 

69 

70 

71@skip_no_h5py 

72def test_serialize_direct_calls_serializer_with_nested_archive(tmp_path: Path) -> None: 

73 """Verify serialize_direct invokes the serializer and returns its 

74 result. 

75 """ 

76 path = str(tmp_path / "test.sdf") 

77 with h5py.File(path, "w") as f: 

78 arch = NdfOutputArchive(f) 

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

80 assert tree.name == "hello" 

81 

82 

83@skip_no_h5py 

84def test_constructor_marks_root_as_ndf(tmp_path: Path) -> None: 

85 """Verify the NdfOutputArchive constructor sets CLASS=NDF on the root 

86 group. 

87 """ 

88 path = str(tmp_path / "test.sdf") 

89 with h5py.File(path, "w") as f: 

90 NdfOutputArchive(f) 

91 with h5py.File(path, "r") as f: 

92 assert f["/"].attrs["CLASS"] == b"NDF" 

93 

94 

95@skip_no_h5py 

96def test_top_level_image_routes_to_data_array(tmp_path: Path) -> None: 

97 """Verify add_array routes a top-level image array to /DATA_ARRAY/DATA.""" 

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

99 path = str(tmp_path / "test.sdf") 

100 with h5py.File(path, "w") as f: 

101 arch = NdfOutputArchive(f) 

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

103 assert ref.source == "ndf:/DATA_ARRAY/DATA" 

104 with h5py.File(path, "r") as f: 

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

106 assert ds.dtype == np.float32 

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

108 assert f["/DATA_ARRAY"].attrs["CLASS"] == b"ARRAY" 

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

110 assert origin.dtype == np.int64 

111 assert origin.shape == (2,) 

112 

113 

114@skip_no_h5py 

115def test_top_level_variance_routes_to_variance(tmp_path: Path) -> None: 

116 """Verify add_array routes a top-level variance array to /VARIANCE/DATA.""" 

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

118 path = str(tmp_path / "test.sdf") 

119 with h5py.File(path, "w") as f: 

120 arch = NdfOutputArchive(f) 

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

122 assert ref.source == "ndf:/VARIANCE/DATA" 

123 with h5py.File(path, "r") as f: 

124 assert f["/VARIANCE"].attrs["CLASS"] == b"ARRAY" 

125 assert f["/VARIANCE/DATA"].dtype == np.float64 

126 

127 

128@skip_no_h5py 

129def test_top_level_compatible_mask_routes_to_quality(tmp_path: Path) -> None: 

130 """Verify add_array routes a 2D uint8 mask to /QUALITY/QUALITY/DATA.""" 

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

132 path = str(tmp_path / "test.sdf") 

133 with h5py.File(path, "w") as f: 

134 arch = NdfOutputArchive(f) 

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

136 assert ref.source == "ndf:/QUALITY/QUALITY/DATA" 

137 with h5py.File(path, "r") as f: 

138 assert f["/QUALITY"].attrs["CLASS"] == b"QUALITY" 

139 assert f["/QUALITY/QUALITY"].attrs["CLASS"] == b"ARRAY" 

140 assert f["/QUALITY/QUALITY/DATA"].dtype == np.uint8 

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

142 assert f["/QUALITY/QUALITY/ORIGIN"].dtype == np.int32 

143 assert f["/QUALITY/QUALITY/ORIGIN"].shape == (2,) 

144 assert f["/QUALITY/QUALITY/BAD_PIXEL"].id.get_type().get_class() == h5py.h5t.BITFIELD 

145 assert not _hds.read_array(f["/QUALITY/QUALITY/BAD_PIXEL"]) 

146 assert f["/QUALITY/BADBITS"][()] == 255 

147 

148 

149@skip_no_h5py 

150def test_top_level_incompatible_mask_routes_to_more_lsst(tmp_path: Path) -> None: 

151 """Verify add_array hoists a 3D uint8 mask to /MORE/LSST/MASK as a sub- 

152 NDF. 

153 """ 

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

155 data[0, 1, 2] = 4 

156 data[1, 2, 3] = 8 

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

158 path = str(tmp_path / "test.sdf") 

159 with h5py.File(path, "w") as f: 

160 arch = NdfOutputArchive(f) 

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

162 assert ref.source == "ndf:/MORE/LSST/MASK/DATA_ARRAY/DATA" 

163 with h5py.File(path, "r") as f: 

164 assert f["/MORE/LSST/MASK"].attrs["CLASS"] == b"NDF" 

165 assert f["/MORE/LSST/MASK/DATA_ARRAY"].attrs["CLASS"] == b"ARRAY" 

166 assert f["/MORE/LSST/MASK/DATA_ARRAY/DATA"].shape == data.shape 

167 assert f["/QUALITY/QUALITY"].attrs["CLASS"] == b"ARRAY" 

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

169 assert f["/QUALITY/BADBITS"][()] == 255 

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

171 assert origin.dtype == np.int64 

172 assert origin.shape == (3,) 

173 

174 

175@skip_no_h5py 

176def test_long_hoisted_component_is_shrunk(tmp_path: Path) -> None: 

177 """Verify HDS component names exceeding DAT__SZNAM are shrunk in the stored 

178 path. 

179 """ 

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

181 path = str(tmp_path / "test.sdf") 

182 with h5py.File(path, "w") as f: 

183 arch = NdfOutputArchive(f) 

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

185 assert ref.source.startswith("ndf:/MORE/LSST/") 

186 assert ref.source.endswith("/DATA_ARRAY/DATA") 

187 with h5py.File(path, "r") as f: 

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

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

190 assert len(component) <= DAT__SZNAM 

191 assert hdf5_path in f 

192 

193 

194@skip_no_h5py 

195def test_long_name_round_trips_through_input_archive(tmp_path: Path) -> None: 

196 """Verify a long-named array can be read back via NdfInputArchive.""" 

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

198 path = str(tmp_path / "test.sdf") 

199 with h5py.File(path, "w") as f: 

200 arch = NdfOutputArchive(f) 

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

202 with NdfInputArchive.open(path) as inp: 

203 read_back = inp.get_array(ref) 

204 np.testing.assert_array_equal(read_back, data) 

205 

206 

207@skip_no_h5py 

208def test_repeated_long_name_gets_distinct_versioned_paths(tmp_path: Path) -> None: 

209 """Verify two identically-named long arrays receive distinct versioned 

210 paths. 

211 """ 

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

213 path = str(tmp_path / "test.sdf") 

214 with h5py.File(path, "w") as f: 

215 arch = NdfOutputArchive(f) 

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

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

218 assert first.source != second.source 

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

220 assert second_leaf.endswith("_2") 

221 with h5py.File(path, "r") as f: 

222 assert first.source[len("ndf:") :] in f 

223 assert second.source[len("ndf:") :] in f 

224 

225 

226@skip_no_h5py 

227def test_nested_array_hoists_as_sub_ndf(tmp_path: Path) -> None: 

228 """Verify nested array names produce a CLASS=NDF sub-structure under 

229 /MORE/LSST. 

230 """ 

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

232 path = str(tmp_path / "test.sdf") 

233 with h5py.File(path, "w") as f: 

234 arch = NdfOutputArchive(f) 

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

236 assert ref.source == "ndf:/MORE/LSST/PSF/COEFFICIENTS/DATA_ARRAY/DATA" 

237 with h5py.File(path, "r") as f: 

238 assert "MORE" in f 

239 assert "LSST" in f["/MORE"] 

240 assert "PSF" in f["/MORE/LSST"] 

241 assert "COEFFICIENTS" in f["/MORE/LSST/PSF"] 

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

243 assert sub.attrs["CLASS"] == b"NDF" 

244 assert sub["DATA_ARRAY"].attrs["CLASS"] == b"ARRAY" 

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

246 origin = sub["DATA_ARRAY/ORIGIN"] 

247 assert origin.dtype == np.int64 

248 assert origin.shape == (data.ndim,) 

249 

250 

251@skip_no_h5py 

252def test_colliding_shrunk_names_raise(tmp_path: Path) -> None: 

253 """Verify add_array raises ValueError when two long names shrink to the 

254 same token. 

255 """ 

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

257 path = str(tmp_path / "test.sdf") 

258 with h5py.File(path, "w") as f: 

259 arch = NdfOutputArchive(f) 

260 with mock.patch.object( 

261 arch._name_shrinker, 

262 "shrink", 

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

264 ): 

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

266 with pytest.raises(ValueError, match="name collision"): 

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

268 

269 

270@skip_no_h5py 

271def test_serialize_pointer_writes_subtree_and_returns_pointer(tmp_path: Path) -> None: 

272 """Verify serialize_pointer stores the sub-tree JSON and returns the 

273 correct pointer. 

274 """ 

275 path = str(tmp_path / "test.sdf") 

276 with h5py.File(path, "w") as f: 

277 arch = NdfOutputArchive(f) 

278 ptr = arch.serialize_pointer( 

279 "psf", 

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

281 key=("psf", 1), 

282 ) 

283 assert ptr.path == "/MORE/LSST/PSF/JSON" 

284 with h5py.File(path, "r") as f: 

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

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

287 assert '"name":"gaussian"' in joined.replace(" ", "") 

288 

289 

290@skip_no_h5py 

291def test_serialize_pointer_caches_by_key(tmp_path: Path) -> None: 

292 """Verify serialize_pointer returns the cached pointer and does not re-run 

293 the serializer. 

294 """ 

295 path = str(tmp_path / "test.sdf") 

296 with h5py.File(path, "w") as f: 

297 arch = NdfOutputArchive(f) 

298 ptr1 = arch.serialize_pointer( 

299 "psf", 

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

301 key=("psf", 1), 

302 ) 

303 ptr2 = arch.serialize_pointer( 

304 "psf", 

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

306 key=("psf", 1), 

307 ) 

308 assert ptr1 == ptr2 

309 with h5py.File(path, "r") as f: 

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

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

312 assert "first" in joined 

313 assert "second" not in joined 

314 

315 

316@skip_no_h5py 

317def test_serialize_pointer_preserves_nested_arrays(tmp_path: Path) -> None: 

318 """Verify serialize_pointer does not clobber nested arrays written by the 

319 serializer. 

320 """ 

321 

322 class TreeWithArray(pydantic.BaseModel): 

323 name: str 

324 data: ArrayReferenceModel 

325 

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

327 path = str(tmp_path / "test.sdf") 

328 with h5py.File(path, "w") as f: 

329 arch = NdfOutputArchive(f) 

330 ptr = arch.serialize_pointer( 

331 "psf", 

332 lambda nested: TreeWithArray( 

333 name="gaussian", 

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

335 ), 

336 key=("psf", 1), 

337 ) 

338 assert ptr.path == "/MORE/LSST/PSF/JSON" 

339 with h5py.File(path, "r") as f: 

340 assert "/MORE/LSST/PSF/JSON" in f 

341 assert "/MORE/LSST/PSF/PARAMETERS/DATA_ARRAY/DATA" in f 

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

343 assert f["/MORE/LSST/PSF"].attrs["CLASS"] == b"PSF" 

344 assert f["/MORE/LSST"].attrs["CLASS"] == b"LSST" 

345 

346 

347@skip_no_h5py 

348def test_serialize_frame_set_records_for_iter(tmp_path: Path) -> None: 

349 """Verify serialize_frame_set records the (FrameSet, pointer) pair for 

350 iter_frame_sets. 

351 """ 

352 frame_set = TinyFrameSet() 

353 path = str(tmp_path / "test.sdf") 

354 with h5py.File(path, "w") as f: 

355 arch = NdfOutputArchive(f) 

356 ptr = arch.serialize_frame_set( 

357 "wcs/pixel_to_sky", 

358 frame_set, 

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

360 key=("frame_set", 1), 

361 ) 

362 assert ptr.path == "/MORE/LSST/WCS/PIXEL_TO_SKY/JSON" 

363 recorded = list(arch.iter_frame_sets()) 

364 assert len(recorded) == 1 

365 assert recorded[0][0] is frame_set 

366 assert recorded[0][1].path == "/MORE/LSST/WCS/PIXEL_TO_SKY/JSON" 

367 

368 

369@skip_no_h5py 

370def test_add_table_returns_inline_table_model(tmp_path: Path) -> None: 

371 """Verify add_table returns an inline table model for a simple astropy 

372 Table. 

373 """ 

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

375 path = str(tmp_path / "test.sdf") 

376 with h5py.File(path, "w") as f: 

377 arch = NdfOutputArchive(f) 

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

379 assert len(model.columns) == 2 

380 assert isinstance(model.columns[0].data, InlineArrayModel) 

381 

382 

383@skip_no_h5py 

384def test_add_structured_array_writes_column_ndfs_with_units(tmp_path: Path) -> None: 

385 """Verify add_structured_array stores each column as a sub-NDF with correct 

386 units. 

387 """ 

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

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

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

391 path = str(tmp_path / "test.sdf") 

392 with h5py.File(path, "w") as f: 

393 arch = NdfOutputArchive(f) 

394 model = arch.add_structured_array( 

395 rec, 

396 name="rec", 

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

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

399 ) 

400 assert len(model.columns) == 2 

401 assert isinstance(model.columns[0].data, ArrayReferenceModel) 

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

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

404 assert col_x.unit == u.m 

405 assert col_y.description == "the y values" 

406 assert col_x.data.source == "ndf:/MORE/LSST/REC/X/DATA_ARRAY/DATA" 

407 assert col_y.data.source == "ndf:/MORE/LSST/REC/Y/DATA_ARRAY/DATA" 

408 with h5py.File(path, "r") as f: 

409 assert f["/MORE/LSST/REC/X"].attrs["CLASS"] == b"NDF" 

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

411 assert f["/MORE/LSST/REC/Y"].attrs["CLASS"] == b"NDF" 

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

413 with NdfInputArchive.open(path) as archive: 

414 recovered = archive.get_structured_array(model) 

415 np.testing.assert_array_equal(recovered, rec) 

416 

417 

418@skip_no_h5py 

419def test_add_single_column_structured_array_uses_table_name(tmp_path: Path) -> None: 

420 """Verify a single-column structured array uses the table path as its NDF 

421 component name. 

422 """ 

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

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

425 path = str(tmp_path / "test.sdf") 

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

427 arch = NdfOutputArchive(f) 

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

429 assert len(model.columns) == 1 

430 column = model.columns[0] 

431 assert isinstance(column.data, ArrayReferenceModel) 

432 assert column.data.source == "ndf:/MORE/LSST/PSF/PIFF/INTERP/SOLUTION/DATA_ARRAY/DATA" 

433 assert column.data.shape == [4] 

434 with h5py.File(path, "r") as f: 

435 assert "PSF" in f["/MORE/LSST"] 

436 assert "PIFF" in f["/MORE/LSST/PSF"] 

437 assert "INTERP" in f["/MORE/LSST/PSF/PIFF"] 

438 assert "SOLUTION" in f["/MORE/LSST/PSF/PIFF/INTERP"] 

439 np.testing.assert_array_equal( 

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

441 rec["solution"], 

442 ) 

443 

444 

445@skip_no_h5py 

446def test_structured_array_long_name_is_shrunk_and_versioned(tmp_path: Path) -> None: 

447 """Verify long structured-array names are shrunk and repeated names get 

448 versioned paths. 

449 """ 

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

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

452 path = str(tmp_path / "test.sdf") 

453 with h5py.File(path, "w") as f: 

454 arch = NdfOutputArchive(f) 

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

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

457 for model in (first, second): 

458 for column in model.columns: 

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

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

461 assert len(component) <= DAT__SZNAM 

462 assert first.columns[0].data.source != second.columns[0].data.source 

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

464 assert second_parent.endswith("_2") 

465 with h5py.File(path, "r") as f: 

466 for model in (first, second): 

467 for column in model.columns: 

468 assert column.data.source[len("ndf:") :] in f 

469 

470 

471@skip_no_h5py 

472def test_write_with_projection_creates_wcs_component(tmp_path: Path) -> None: 

473 """Verify write() creates a /WCS/DATA component when the image has a 

474 sky_projection. 

475 """ 

476 rng = np.random.default_rng(42) 

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

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

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

480 image = Image( 

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

482 bbox=bbox, 

483 sky_projection=sky_projection, 

484 ) 

485 path = str(tmp_path / "test.sdf") 

486 write(image, path) 

487 with h5py.File(path, "r") as f: 

488 assert "WCS" in f 

489 assert f["/WCS"].attrs["CLASS"] == b"WCS" 

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

491 assert wcs_data.dtype == np.dtype("|S32") 

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

493 assert all(record[0] in {" ", "+"} for record in records) 

494 assert not any(record.startswith("#") for record in records) 

495 text = _hds.decode_ndf_ast_data(records) 

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

497 assert any(s.startswith("Begin FrameSet") for s in stripped) 

498 assert any(s.startswith("End FrameSet") for s in stripped) 

499 assert 'Domain = "GRID"' in stripped 

500 assert 'Domain = "PIXEL"' in stripped 

501 assert "Sft1 = -19" in stripped 

502 assert "Sft2 = -9" in stripped 

503 

504 

505@skip_no_h5py 

506def test_write_without_projection_omits_wcs_component(tmp_path: Path) -> None: 

507 """Verify write() omits /WCS when the image has no sky_projection.""" 

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

509 path = str(tmp_path / "test.sdf") 

510 write(image, path) 

511 with h5py.File(path, "r") as f: 

512 assert "WCS" not in f 

513 

514 

515@skip_no_h5py 

516def test_mask_sub_ndf_gets_3d_wcs(tmp_path: Path) -> None: 

517 """Verify an incompatible mask hoisted to /MORE/LSST/MASK carries a 3D 

518 /WCS. 

519 """ 

520 rng = np.random.default_rng(42) 

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

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

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

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

525 image = Image( 

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

527 bbox=bbox, 

528 sky_projection=sky_projection, 

529 ) 

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

531 path = str(tmp_path / "test.sdf") 

532 write(masked, path) 

533 with h5py.File(path, "r") as f: 

534 assert "WCS" in f 

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

536 assert "MASK" in f["/MORE/LSST"] 

537 assert "WCS" in f["/MORE/LSST/MASK"] 

538 assert f["/MORE/LSST/MASK/WCS"].attrs["CLASS"] == b"WCS" 

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

540 assert top_lines != mask_lines 

541 mask_text = _hds.decode_ndf_ast_data(mask_lines) 

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

543 assert "Naxes = 3" in stripped 

544 assert 'Domain = "GRID"' in stripped 

545 assert 'Domain = "PIXEL"' in stripped 

546 assert "Sft1 = -19" in stripped 

547 assert "Sft2 = -9" in stripped 

548 assert "Sft3 = 1" in stripped 

549 assert "Begin CmpFrame" in stripped 

550 assert "Begin SkyFrame" in stripped 

551 assert 'Domain = "MASK"' in stripped 

552 assert "Begin CmpMap" in stripped 

553 assert "Series = 0" in stripped 

554 

555 

556@skip_no_h5py 

557def test_mask_sub_ndf_no_wcs_when_image_has_no_projection(tmp_path: Path) -> None: 

558 """Verify /MORE/LSST/MASK does not carry /WCS when the image has no 

559 sky_projection. 

560 """ 

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

562 masked = MaskedImage( 

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

564 mask_schema=MaskSchema(planes), 

565 ) 

566 path = str(tmp_path / "test.sdf") 

567 write(masked, path) 

568 with h5py.File(path, "r") as f: 

569 assert "WCS" not in f 

570 assert "MASK" in f["/MORE/LSST"] 

571 assert "WCS" not in f["/MORE/LSST/MASK"] 

572 

573 

574@skip_no_h5py 

575def test_write_image_produces_valid_layout(tmp_path: Path) -> None: 

576 """Verify write() produces a valid NDF layout with DATA_ARRAY, ORIGIN, and 

577 /MORE/LSST/JSON. 

578 """ 

579 image = Image( 

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

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

582 ) 

583 path = str(tmp_path / "test.sdf") 

584 tree = write(image, path) 

585 assert tree is not None 

586 with h5py.File(path, "r") as f: 

587 assert f["/"].attrs["CLASS"] == b"NDF" 

588 assert "HDS_ROOT_NAME" in f["/"].attrs 

589 assert f["/DATA_ARRAY"].attrs["CLASS"] == b"ARRAY" 

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

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

592 assert origin.dtype == np.int64 

593 assert len(origin) == 2 

594 assert not (origin == 0).all() 

595 assert "MORE" in f 

596 assert "LSST" in f["/MORE"] 

597 assert "JSON" in f["/MORE/LSST"] 

598 

599 

600@skip_no_h5py 

601def test_write_image_preserves_opaque_fits_metadata(tmp_path: Path) -> None: 

602 """Verify write() stores opaque FITS headers in /MORE/FITS and they survive 

603 a round-trip. 

604 """ 

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

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

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

608 long_value = "x" * 100 

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

610 opaque = FitsOpaqueMetadata() 

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

612 image._opaque_metadata = opaque 

613 path = str(tmp_path / "test.sdf") 

614 write(image, path) 

615 with h5py.File(path, "r") as f: 

616 assert "FITS" in f["/MORE"] 

617 cards = [c.decode("ascii").rstrip(" ") for c in f["/MORE/FITS"][()]] 

618 assert any(c.startswith("FOO") for c in cards) 

619 assert any(c.startswith("CONTINUE") for c in cards) 

620 assert all(len(c.encode("ascii")) <= 80 for c in cards) 

621 result = read_archive(path, Image) 

622 recovered = result._opaque_metadata.headers[ExtensionKey()] 

623 assert recovered["LONGSTR"] == long_value 

624 

625 

626@skip_no_h5py 

627def test_write_image_main_json_round_trips_back(tmp_path: Path) -> None: 

628 """Verify the main JSON tree at /MORE/LSST/JSON matches write()'s returned 

629 ArchiveTree. 

630 """ 

631 image = Image(np.arange(6, dtype=np.float32).reshape(2, 3)) 

632 path = str(tmp_path / "test.sdf") 

633 tree = write(image, path) 

634 with h5py.File(path, "r") as f: 

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

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

637 recovered = json.loads(joined) 

638 assert json.loads(tree.model_dump_json()) == recovered 

639 

640 

641@skip_no_h5py 

642def test_write_image_with_unit_creates_units_component(tmp_path: Path) -> None: 

643 """Verify write() creates a /UNITS component and it round-trips back as the 

644 correct unit. 

645 """ 

646 image = Image(np.arange(6, dtype=np.float32).reshape(2, 3), unit=u.ct) 

647 path = str(tmp_path / "test.sdf") 

648 write(image, path) 

649 with h5py.File(path, "r") as f: 

650 assert "UNITS" in f 

651 assert f["/UNITS"].shape == () 

652 assert f["/UNITS"][()].decode("ascii").rstrip(" ") == "count" 

653 result = read_archive(path, Image) 

654 assert result.unit == u.ct