Coverage for tests/test_fields.py: 92%

318 statements  

« prev     ^ index     » next       coverage.py v7.15.2, created at 2026-07-16 15:24 -0700

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 os 

15from collections.abc import Callable 

16from typing import Any 

17 

18import astropy.units as u 

19import numpy as np 

20import pytest 

21 

22from lsst.images import XY, YX, Box, Image 

23from lsst.images.fields import ( 

24 BaseField, 

25 ChebyshevField, 

26 ProductField, 

27 SplineField, 

28 SumField, 

29 field_from_legacy, 

30 field_from_legacy_background, 

31) 

32from lsst.images.tests import ( 

33 RoundtripFits, 

34 assert_close, 

35 assert_images_equal, 

36 compare_field_to_legacy, 

37) 

38 

39try: 

40 from lsst.afw.image import MaskedImageF 

41 from lsst.afw.math import BackgroundList as LegacyBackgroundList 

42 from lsst.afw.math import BackgroundMI 

43 from lsst.afw.math import Chebyshev1Function2D as LegacyChebyshev1Function2D 

44 from lsst.afw.math import ChebyshevBoundedField as LegacyChebyshevBoundedField 

45 from lsst.afw.math import ProductBoundedField as LegacyProductBoundedField 

46 from lsst.geom import Box2D as LegacyBox2D 

47 

48 HAVE_LEGACY = True 

49except ImportError: 

50 HAVE_LEGACY = False 

51 type LegacyBackgroundList = Any # type: ignore[no-redef] 

52 type LegacyBox2D = Any # type: ignore[no-redef] 

53 type LegacyChebyshev1Function2D = Any # type: ignore[no-redef] 

54 type LegacyChebyshevBoundedField = Any # type: ignore[no-redef] 

55 type LegacyProductBoundedField = Any # type: ignore[no-redef] 

56 

57 

58EXTERNAL_DATA_DIR = os.environ.get("TESTDATA_IMAGES_DIR", None) 

59TEST_BOX = Box.factory[6:32, -7:26] 

60 

61 

62@pytest.fixture(scope="session") 

63def legacy_visit_background() -> LegacyBackgroundList: 

64 """Load and return an `lsst.afw.math.BackgroundList`. 

65 

66 Skips if TESTDATA_IMAGES_DIR is unset or lsst.afw.math is unavailable. 

67 """ 

68 if EXTERNAL_DATA_DIR is None: 68 ↛ 70line 68 didn't jump to line 70 because the condition on line 68 was always true

69 pytest.skip("TESTDATA_IMAGES_DIR is not in the environment.") 

70 if not HAVE_LEGACY: 

71 pytest.skip("This test requires lsst.afw.math to be importable.") 

72 filename = os.path.join(EXTERNAL_DATA_DIR, "dp2", "legacy", "visit_image_background.fits") 

73 return LegacyBackgroundList.readFits(filename) 

74 

75 

76def _make_test_chebyshev_field() -> ChebyshevField: 

77 return ChebyshevField(TEST_BOX, np.array([[0.5, -0.25], [0.40, 0.0]])) 

78 

79 

80def _make_test_spline_field() -> SplineField: 

81 rng = np.random.default_rng(10) 

82 return SplineField( 

83 TEST_BOX, 

84 rng.standard_normal(size=(6, 7)), 

85 y=TEST_BOX.y.linspace(6), 

86 x=TEST_BOX.x.linspace(7), 

87 ) 

88 

89 

90def _make_test_product_field() -> ProductField: 

91 return ProductField([_make_test_spline_field(), _make_test_chebyshev_field()]) 

92 

93 

94def _make_test_sum_field() -> SumField: 

95 return SumField([_make_test_spline_field(), _make_test_chebyshev_field()]) 

96 

97 

98def _make_legacy_chebyshev_field() -> LegacyChebyshevBoundedField: 

99 if not HAVE_LEGACY: 99 ↛ 100line 99 didn't jump to line 100 because the condition on line 99 was never true

100 pytest.skip("Legacy LSST packages could not be imported.") 

101 rng = np.random.default_rng(10) 

102 cheby_coeffs = rng.random((6, 3)) 

103 return LegacyChebyshevBoundedField(TEST_BOX.to_legacy(), cheby_coeffs) 

104 

105 

106def _make_legacy_product_field() -> LegacyProductBoundedField: 

107 if not HAVE_LEGACY: 107 ↛ 108line 107 didn't jump to line 108 because the condition on line 107 was never true

108 pytest.skip("Legacy LSST packages could not be imported.") 

109 rng = np.random.default_rng(11) 

110 cheby2 = LegacyChebyshevBoundedField(TEST_BOX.to_legacy(), rng.standard_normal(size=(2, 2))) 

111 return LegacyProductBoundedField([_make_legacy_chebyshev_field(), cheby2]) 

112 

113 

114def check_evaluation_consistency(field: BaseField) -> None: 

115 """Check that __call__ and render agree.""" 

116 image_1 = field.render() 

117 p = field.bounds.bbox.meshgrid() 

118 image_2 = Image(field(x=p.x, y=p.y), bbox=field.bounds.bbox, unit=field.unit) 

119 assert_images_equal(image_1, image_2) 

120 scaled_field = field * 2.0 

121 image_3 = scaled_field.render() 

122 image_3.array *= 0.5 

123 assert_images_equal(image_1, image_3) 

124 image_4 = field.render(Box.factory[9:11, -3:1]) 

125 assert_images_equal(image_4, image_1[image_4.bbox]) 

126 

127 

128@pytest.mark.parametrize( 

129 "factory", 

130 [ 

131 _make_test_chebyshev_field, 

132 _make_test_spline_field, 

133 _make_test_product_field, 

134 _make_test_sum_field, 

135 ], 

136) 

137def test_evaluation_consistency(factory: Callable[[], BaseField]) -> None: 

138 """Test that __call__ and render agree.""" 

139 check_evaluation_consistency(factory()) 

140 

141 

142@pytest.mark.parametrize( 

143 "factory", 

144 [ 

145 _make_test_chebyshev_field, 

146 _make_test_spline_field, 

147 _make_test_product_field, 

148 _make_test_sum_field, 

149 ], 

150) 

151def test_units(factory: Callable[[], BaseField]) -> None: 

152 """Test that fields correctly propagats and check units.""" 

153 field = factory() 

154 assert field.unit is None 

155 with_units_1 = field * u.nJy 

156 assert with_units_1.unit == u.nJy 

157 assert with_units_1(x=np.array([0.0]), y=np.array([10.0]), quantity=True).unit == u.nJy 

158 image_1 = with_units_1.render(bbox=Box.factory[10:12, 0:3]) 

159 assert image_1.unit == u.nJy 

160 assert (with_units_1 * 2.0).unit == u.nJy 

161 assert (with_units_1 / u.arcsec**2).unit == u.nJy / u.arcsec**2 

162 

163 

164def test_chebyshev_call_limits() -> None: 

165 """Test that ChebyshevField evaluates correctly at first order at the 

166 corners of its box. 

167 """ 

168 cheby = _make_test_chebyshev_field() 

169 result = cheby(x=np.array([-7.5, 25.5, 25.5, -7.5]), y=np.array([5.5, 5.5, 31.5, 31.5])) 

170 assert result[0] == 0.5 + 0.25 - 0.4 # [x=-1, y=-1] after remap 

171 assert result[1] == 0.5 - 0.25 - 0.4 # [x=1, y=-1] after remap 

172 assert result[2] == 0.5 - 0.25 + 0.4 # [x=1, y=1] after remap 

173 assert result[3] == 0.5 + 0.25 + 0.4 # [x=-1, y=1] after remap 

174 

175 

176def test_chebyshev_attributes() -> None: 

177 """Test the basic properties of a ChebyshevField.""" 

178 cheby = _make_test_chebyshev_field() 

179 assert cheby.bounds == TEST_BOX 

180 assert cheby.unit is None 

181 assert cheby.x_order == 1 

182 assert cheby.y_order == 1 

183 assert cheby.order == 1 

184 np.testing.assert_array_equal(cheby.coefficients, np.array([[0.5, -0.25], [0.40, 0.0]])) 

185 

186 

187def test_chebyshev_fit() -> None: 

188 """Test that ChebyshevField.fit recovers the original coefficients with 

189 zero residuals. 

190 """ 

191 cheby = _make_test_chebyshev_field() 

192 rng = np.random.default_rng(22) 

193 data_image = cheby.render() 

194 cheby2 = ChebyshevField.fit(TEST_BOX, data_image.array, order=1, y=TEST_BOX.y.arange, x=TEST_BOX.x.arange) 

195 assert cheby2.bounds == TEST_BOX 

196 np.testing.assert_array_almost_equal(cheby2.coefficients, cheby.coefficients) 

197 # Fit to order 2 in x (will get us extra zero-valued coefficients): 

198 cheby3 = ChebyshevField.fit( 

199 TEST_BOX, data_image.array, x_order=2, y_order=1, y=TEST_BOX.y.arange, x=TEST_BOX.x.arange 

200 ) 

201 assert cheby3.bounds == TEST_BOX 

202 np.testing.assert_array_almost_equal( 

203 cheby3.coefficients, 

204 np.array([[0.5, -0.25, 0.0], [0.40, 0.0, 0.0]], dtype=np.float64), 

205 ) 

206 # Fit with triangular=False: 

207 cheby4 = ChebyshevField.fit( 

208 TEST_BOX, data_image.array, order=1, y=TEST_BOX.y.arange, x=TEST_BOX.x.arange, triangular=False 

209 ) 

210 assert cheby4.bounds == TEST_BOX 

211 np.testing.assert_array_almost_equal(cheby4.coefficients, cheby.coefficients) 

212 # Fit with weights. 

213 cheby5 = ChebyshevField.fit( 

214 TEST_BOX, 

215 data_image.array, 

216 order=1, 

217 y=TEST_BOX.y.arange, 

218 x=TEST_BOX.x.arange, 

219 weight=rng.uniform(0.8, 1.2, size=data_image.array.shape), 

220 ) 

221 assert cheby5.bounds == TEST_BOX 

222 np.testing.assert_array_almost_equal(cheby5.coefficients, cheby.coefficients) 

223 # Fit to a Quantity. 

224 cheby6 = ChebyshevField.fit( 

225 TEST_BOX, data_image.array * u.nJy, order=1, y=TEST_BOX.y.arange, x=TEST_BOX.x.arange 

226 ) 

227 assert cheby6.bounds == TEST_BOX 

228 assert cheby6.unit == u.nJy 

229 np.testing.assert_array_almost_equal(cheby6.coefficients, cheby.coefficients) 

230 # Fit with units provided separately. 

231 cheby7 = ChebyshevField.fit( 

232 TEST_BOX, data_image.array, order=1, y=TEST_BOX.y.arange, x=TEST_BOX.x.arange, unit=u.nJy 

233 ) 

234 assert cheby7.bounds == TEST_BOX 

235 assert cheby7.unit == u.nJy 

236 np.testing.assert_array_almost_equal(cheby7.coefficients, cheby.coefficients) 

237 # Fit with x and y labeling every data point. 

238 m = TEST_BOX.meshgrid() 

239 cheby8 = ChebyshevField.fit(TEST_BOX, data_image.array, order=1, y=m.y, x=m.x) 

240 assert cheby8.bounds == TEST_BOX 

241 np.testing.assert_array_almost_equal(cheby8.coefficients, cheby.coefficients) 

242 # Fit with x and y labeling every data point plus weights. 

243 cheby9 = ChebyshevField.fit( 

244 TEST_BOX, 

245 data_image.array, 

246 order=1, 

247 y=m.y, 

248 x=m.x, 

249 weight=rng.uniform(0.8, 1.2, size=data_image.array.shape), 

250 ) 

251 assert cheby9.bounds == TEST_BOX 

252 np.testing.assert_array_almost_equal(cheby9.coefficients, cheby.coefficients) 

253 # Fit with one data point replaced by NaN (should be ignored). 

254 new_data = data_image.array.copy() 

255 new_data[5, 7] = np.nan 

256 cheby10 = ChebyshevField.fit(TEST_BOX, new_data, order=1, y=TEST_BOX.y.arange, x=TEST_BOX.x.arange) 

257 assert cheby10.bounds == TEST_BOX 

258 np.testing.assert_array_almost_equal(cheby10.coefficients, cheby.coefficients) 

259 

260 

261def test_spline_knot_evaluation() -> None: 

262 """Test that SplineField evaluates to its input data at the knot 

263 positions. 

264 """ 

265 spline = _make_test_spline_field() 

266 xv, yv = np.meshgrid(spline.x, spline.y) 

267 result = spline(x=xv, y=yv) 

268 np.testing.assert_array_almost_equal(result, spline.data) 

269 

270 

271def test_product_evaluation() -> None: 

272 """Test that ProductField.__call__ against direct calls to its operands.""" 

273 product = _make_test_product_field() 

274 cheby = _make_test_chebyshev_field() 

275 spline = _make_test_spline_field() 

276 xv, yv = TEST_BOX.meshgrid(n=3) 

277 z = product(x=xv, y=yv) 

278 np.testing.assert_array_equal(z, cheby(x=xv, y=yv) * spline(x=xv, y=yv)) 

279 

280 

281def test_product_units() -> None: 

282 """Test that ProductField correctly propagates and combines units.""" 

283 cheby = _make_test_chebyshev_field() 

284 spline = _make_test_spline_field() 

285 assert ProductField([cheby, spline * u.nJy]).unit == u.nJy 

286 assert ProductField([cheby * u.nJy, spline]).unit == u.nJy 

287 assert ProductField([cheby * u.nJy, spline / u.arcsec**2]).unit == u.nJy / u.arcsec**2 

288 

289 

290def test_sum_evaluation() -> None: 

291 """Test that SumField.__call__ against direct calls to its operands.""" 

292 sum_ = _make_test_sum_field() 

293 cheby = _make_test_chebyshev_field() 

294 spline = _make_test_spline_field() 

295 xv, yv = TEST_BOX.meshgrid(n=3) 

296 z = sum_(x=xv, y=yv) 

297 np.testing.assert_array_equal(z, cheby(x=xv, y=yv) + spline(x=xv, y=yv)) 

298 

299 

300def test_sum_units() -> None: 

301 """Test that SumField correctly propagates units and raises on incompatible 

302 units. 

303 """ 

304 cheby = _make_test_chebyshev_field() 

305 spline = _make_test_spline_field() 

306 with pytest.raises(u.UnitConversionError): 

307 SumField([cheby, spline * u.nJy]) 

308 with pytest.raises(u.UnitConversionError): 

309 SumField([cheby * u.nJy, spline]) 

310 # Test a SumField where the operands have different but compatible units. 

311 mixed = SumField([cheby * u.rad, spline * u.deg]) 

312 assert mixed.unit == u.rad 

313 small_box = Box.factory[10:12, 2:5] 

314 cheby_render = cheby.render(small_box) 

315 spline_render = spline.render(small_box) 

316 mixed_render = mixed.render(small_box) 

317 assert mixed_render.unit == u.rad 

318 np.testing.assert_array_almost_equal( 

319 mixed_render.array, cheby_render.array + (spline_render.array * np.pi / 180.0) 

320 ) 

321 check_evaluation_consistency(mixed) 

322 

323 

324def test_chebyshev_roundtrip() -> None: 

325 """Test converting ChebyshevField from and to legacy, and serializing it 

326 in between. 

327 """ 

328 legacy_cheby = _make_legacy_chebyshev_field() 

329 cheby = field_from_legacy(legacy_cheby) 

330 assert isinstance(cheby, ChebyshevField) 

331 compare_field_to_legacy(cheby, legacy_cheby, subimage_bbox=Box.factory[8:12, -3:2]) 

332 with RoundtripFits(cheby) as roundtrip: 

333 pass 

334 compare_field_to_legacy(roundtrip.result, legacy_cheby, subimage_bbox=Box.factory[8:12, -3:2]) 

335 compare_field_to_legacy(roundtrip.result, cheby.to_legacy(), subimage_bbox=Box.factory[8:12, -3:2]) 

336 

337 

338def test_product_roundtrip() -> None: 

339 """Test converting ProductField from and to legacy, and serializing it 

340 in between. 

341 """ 

342 legacy_product = _make_legacy_product_field() 

343 product = field_from_legacy(legacy_product) 

344 assert isinstance(product, ProductField) 

345 compare_field_to_legacy(product, legacy_product, subimage_bbox=Box.factory[8:12, -3:2]) 

346 with RoundtripFits(product) as roundtrip: 

347 pass 

348 compare_field_to_legacy(roundtrip.result, legacy_product, subimage_bbox=Box.factory[8:12, -3:2]) 

349 compare_field_to_legacy(roundtrip.result, product.to_legacy(), subimage_bbox=Box.factory[8:12, -3:2]) 

350 

351 

352def test_spline_simple() -> None: 

353 """Test SplineField against BackgroundMI with no missing data.""" 

354 if not HAVE_LEGACY: 354 ↛ 355line 354 didn't jump to line 355 because the condition on line 354 was never true

355 pytest.skip("Legacy LSST packages could not be imported.") 

356 rng = np.random.default_rng(23) 

357 bins = MaskedImageF(Box.factory[0:5, 0:6].to_legacy()) 

358 bins.image.array[:, :] = rng.standard_normal(bins.image.array.shape) 

359 bins.variance.array[::] = 1.0 

360 legacy_bg = BackgroundMI(TEST_BOX.to_legacy(), bins) 

361 spline = field_from_legacy_background(legacy_bg) 

362 render_bbox = TEST_BOX.padded(-3) 

363 assert_images_equal( 

364 spline.render(render_bbox), 

365 Image.from_legacy( 

366 legacy_bg.getImageF(render_bbox.to_legacy(), legacy_bg.getBackgroundControl().getInterpStyle()) 

367 ), 

368 rtol=1e-7, 

369 ) 

370 

371 

372def test_spline_one_nan() -> None: 

373 """Test SplineField against BackgroundMI with one missing data point.""" 

374 if not HAVE_LEGACY: 374 ↛ 375line 374 didn't jump to line 375 because the condition on line 374 was never true

375 pytest.skip("Legacy LSST packages could not be imported.") 

376 rng = np.random.default_rng(24) 

377 bins = MaskedImageF(Box.factory[0:7, 0:6].to_legacy()) 

378 bins.image.array[:, :] = rng.standard_normal(bins.image.array.shape) 

379 bins.image.array[3, 2] = np.nan 

380 bins.variance.array[::] = 1.0 

381 legacy_bg = BackgroundMI(TEST_BOX.to_legacy(), bins) 

382 spline = field_from_legacy_background(legacy_bg) 

383 render_bbox = TEST_BOX.padded(-3) 

384 assert_images_equal( 

385 spline.render(render_bbox), 

386 Image.from_legacy( 

387 legacy_bg.getImageF( 

388 render_bbox.to_legacy(), 

389 legacy_bg.getBackgroundControl().getInterpStyle(), 

390 ) 

391 ), 

392 rtol=1e-7, 

393 ) 

394 

395 

396def test_chebyshev1_function2() -> None: 

397 """Verify ChebyshevField.from_legacy_function2 and to_legacy_function2 

398 round-trip. 

399 """ 

400 if not HAVE_LEGACY: 400 ↛ 401line 400 didn't jump to line 401 because the condition on line 400 was never true

401 pytest.skip("Legacy LSST packages could not be imported.") 

402 rng = np.random.default_rng(25) 

403 legacy_func2a = LegacyChebyshev1Function2D(4, LegacyBox2D(TEST_BOX.to_legacy())) 

404 legacy_func2a.setParameters(rng.standard_normal(legacy_func2a.getNParameters())) 

405 field = ChebyshevField.from_legacy_function2(legacy_func2a) 

406 legacy_func2b = field.to_legacy_function2() 

407 assert field.bounds == TEST_BOX 

408 xy_array = TEST_BOX.meshgrid(4) 

409 z_array = field(x=xy_array.x, y=xy_array.y) 

410 for z, x, y in zip(z_array.flat, xy_array.x.flat, xy_array.y.flat): 

411 assert_close(legacy_func2a(x, y), z) 

412 assert_close(legacy_func2b(x, y), z) 

413 

414 

415def test_visit_background(legacy_visit_background: LegacyBackgroundList) -> None: 

416 """Test field_from_legacy_background against a real visit image 

417 background. 

418 """ 

419 bg_field = field_from_legacy_background(legacy_visit_background) 

420 assert_images_equal(bg_field.render(), Image.from_legacy(legacy_visit_background.getImage()), rtol=1e-6) 

421 

422 

423# Scalar x/y inside TEST_BOX used by the broadcasting/scalar tests below. 

424_SCALAR_X = 9.0 

425_SCALAR_Y = 18.5 

426_SCALAR_X_INT = 9 

427_SCALAR_Y_INT = 18 

428 

429 

430@pytest.mark.parametrize( 

431 "factory", 

432 [ 

433 _make_test_chebyshev_field, 

434 _make_test_spline_field, 

435 _make_test_product_field, 

436 _make_test_sum_field, 

437 ], 

438) 

439def test_call_scalar(factory: Callable[[], BaseField]) -> None: 

440 """Verify that __call__ accepts scalar x/y and returns a scalar float, 

441 and returns a 0-d Quantity when quantity=True. 

442 """ 

443 field = factory() 

444 # quantity=False: Python int and float scalars should return float. 

445 result_float = field(x=_SCALAR_X, y=_SCALAR_Y) 

446 assert type(result_float) is float 

447 result_int = field(x=_SCALAR_X_INT, y=_SCALAR_Y_INT) 

448 assert type(result_int) is float 

449 # Values must match the corresponding single-element array call. 

450 ref_float = field(x=np.array([_SCALAR_X]), y=np.array([_SCALAR_Y]))[0] 

451 assert result_float == ref_float 

452 ref_int = field(x=np.array([float(_SCALAR_X_INT)]), y=np.array([float(_SCALAR_Y_INT)]))[0] 

453 assert result_int == ref_int 

454 # quantity=True: scalar input should yield a 0-d Quantity. 

455 field_with_units = field * u.nJy 

456 result_q = field_with_units(x=_SCALAR_X, y=_SCALAR_Y, quantity=True) 

457 assert isinstance(result_q, u.Quantity) 

458 assert result_q.shape == () 

459 assert result_q.unit == u.nJy 

460 ref_q = field_with_units(x=np.array([_SCALAR_X]), y=np.array([_SCALAR_Y]), quantity=True) 

461 assert result_q.value == ref_q[0].value 

462 

463 

464@pytest.mark.parametrize( 

465 "factory", 

466 [ 

467 _make_test_chebyshev_field, 

468 _make_test_spline_field, 

469 _make_test_product_field, 

470 _make_test_sum_field, 

471 ], 

472) 

473def test_call_array_like_and_integer_input(factory: Callable[[], BaseField]) -> None: 

474 """Verify that __call__ accepts Python lists and integer arrays, returning 

475 float64 ndarray results consistent with float64 array input. 

476 """ 

477 field = factory() 

478 xv = [_SCALAR_X, _SCALAR_X + 2.0, _SCALAR_X + 4.0] 

479 yv = [_SCALAR_Y, _SCALAR_Y + 1.0, _SCALAR_Y + 2.0] 

480 # Python list input should return an ndarray. 

481 result_list = field(x=xv, y=yv) 

482 assert isinstance(result_list, np.ndarray) 

483 ref = field(x=np.array(xv), y=np.array(yv)) 

484 np.testing.assert_array_equal(result_list, ref) 

485 # Integer array input should not raise and should return float64. 

486 xi = np.array([_SCALAR_X_INT, _SCALAR_X_INT + 2, _SCALAR_X_INT + 4], dtype=np.int32) 

487 yi = np.array([_SCALAR_Y_INT, _SCALAR_Y_INT + 1, _SCALAR_Y_INT + 2], dtype=np.int32) 

488 result_int = field(x=xi, y=yi) 

489 assert isinstance(result_int, np.ndarray) 

490 assert result_int.dtype == np.float64 

491 ref_int = field(x=xi.astype(np.float64), y=yi.astype(np.float64)) 

492 np.testing.assert_array_equal(result_int, ref_int) 

493 

494 

495@pytest.mark.parametrize( 

496 "factory", 

497 [ 

498 _make_test_chebyshev_field, 

499 _make_test_spline_field, 

500 _make_test_product_field, 

501 _make_test_sum_field, 

502 ], 

503) 

504def test_call_broadcast(factory: Callable[[], BaseField]) -> None: 

505 """Verify that __call__ broadcasts x and y like a NumPy ufunc, in both 

506 1-D and 2-D cases. 

507 """ 

508 field = factory() 

509 xv = np.linspace(_SCALAR_X, _SCALAR_X + 4.0, 5) 

510 yv = np.linspace(_SCALAR_Y, _SCALAR_Y + 3.0, 4) 

511 # 1-D broadcast: array x, scalar y. 

512 result_1d = field(x=xv, y=_SCALAR_Y) 

513 assert isinstance(result_1d, np.ndarray) 

514 assert result_1d.shape == xv.shape 

515 ref_1d = field(x=xv, y=np.full_like(xv, _SCALAR_Y)) 

516 np.testing.assert_array_equal(result_1d, ref_1d) 

517 # 2-D broadcast: column x (M,1), row y (1,N) -> (M,N). 

518 x2d = xv[:, np.newaxis] # shape (5, 1) 

519 y2d = yv[np.newaxis, :] # shape (1, 4) 

520 result_2d = field(x=x2d, y=y2d) 

521 assert isinstance(result_2d, np.ndarray) 

522 assert result_2d.shape == (xv.size, yv.size) 

523 # Values must match the fully expanded meshgrid call. 

524 xmesh, ymesh = np.meshgrid(xv, yv, indexing="ij") 

525 ref_2d = field(x=xmesh, y=ymesh) 

526 np.testing.assert_array_equal(result_2d, ref_2d) 

527 

528 

529@pytest.mark.parametrize( 

530 "factory", 

531 [ 

532 _make_test_chebyshev_field, 

533 _make_test_spline_field, 

534 _make_test_product_field, 

535 _make_test_sum_field, 

536 ], 

537) 

538def test_call_float32_input(factory: Callable[[], BaseField]) -> None: 

539 """Verify that float32 inputs produce float64 output with correct values, 

540 without silent precision loss. 

541 """ 

542 field = factory() 

543 x32 = np.array([_SCALAR_X, _SCALAR_X + 2.0, _SCALAR_X + 4.0], dtype=np.float32) 

544 y32 = np.array([_SCALAR_Y, _SCALAR_Y + 1.0, _SCALAR_Y + 2.0], dtype=np.float32) 

545 result = field(x=x32, y=y32) 

546 assert isinstance(result, np.ndarray) 

547 assert result.dtype == np.float64 

548 ref = field(x=x32.astype(np.float64), y=y32.astype(np.float64)) 

549 np.testing.assert_array_equal(result, ref) 

550 

551 

552@pytest.mark.parametrize( 

553 "factory", 

554 [ 

555 _make_test_chebyshev_field, 

556 _make_test_spline_field, 

557 _make_test_product_field, 

558 _make_test_sum_field, 

559 ], 

560) 

561def test_call_xy_yx(factory: Callable[[], BaseField]) -> None: 

562 """Verify that __call__ accepts XY and YX positional arguments, producing 

563 results identical to the equivalent x=/y= keyword call. 

564 """ 

565 field = factory() 

566 # Scalar XY and YX — return type must be float and values must match. 

567 ref_scalar = field(x=_SCALAR_X, y=_SCALAR_Y) 

568 assert field(XY(x=_SCALAR_X, y=_SCALAR_Y)) == ref_scalar 

569 assert field(YX(y=_SCALAR_Y, x=_SCALAR_X)) == ref_scalar 

570 # Array XY and YX — return type must be ndarray and values must match. 

571 xv = np.array([_SCALAR_X, _SCALAR_X + 2.0, _SCALAR_X + 4.0]) 

572 yv = np.array([_SCALAR_Y, _SCALAR_Y + 1.0, _SCALAR_Y + 2.0]) 

573 ref_array = field(x=xv, y=yv) 

574 np.testing.assert_array_equal(field(XY(xv, yv)), ref_array) 

575 np.testing.assert_array_equal(field(YX(yv, xv)), ref_array) 

576 # quantity=True still works alongside an XY positional argument. 

577 field_with_units = field * u.nJy 

578 ref_q = field_with_units(x=_SCALAR_X, y=_SCALAR_Y, quantity=True) 

579 result_q = field_with_units(XY(x=_SCALAR_X, y=_SCALAR_Y), quantity=True) 

580 assert isinstance(result_q, u.Quantity) 

581 assert result_q.value == ref_q.value 

582 # Mixing a positional point with explicit x= or y= must raise TypeError. 

583 with pytest.raises(TypeError): 

584 field(XY(x=_SCALAR_X, y=_SCALAR_Y), x=_SCALAR_X) 

585 with pytest.raises(TypeError): 

586 field(YX(y=_SCALAR_Y, x=_SCALAR_X), y=_SCALAR_Y)