Coverage for python/lsst/images/fields/_base.py: 35%

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 

14__all__ = ("BaseField",) 

15 

16from abc import ABC, abstractmethod 

17from typing import TYPE_CHECKING, Any, Literal, Self, assert_type, cast, overload 

18 

19import astropy.units 

20import numpy as np 

21import numpy.typing as npt 

22 

23from .._geom import XY, YX, Bounds, Box 

24from .._image import Image 

25from ..describe import DescribableMixin, DescribeOptions, FieldRole, Report, ReportField 

26 

27if TYPE_CHECKING: 

28 try: 

29 from lsst.afw.image import PhotoCalib as LegacyPhotoCalib 

30 from lsst.afw.math import BoundedField as LegacyBoundedField 

31 except ImportError: 

32 type LegacyBoundedField = Any # type: ignore[no-redef] 

33 type LegacyPhotoCalib = Any # type: ignore[no-redef] 

34 

35 

36class BaseField(DescribableMixin, ABC): 

37 """An abstract base class for parametric or interpolated 2-d functions, 

38 generally representing some sort of calculated image. 

39 

40 Notes 

41 ----- 

42 The field hierarchy is closed to the types in this package, so we can 

43 enumerate all of the serializations and avoid any kind of extension system. 

44 All field types are immutable. 

45 

46 Field types implement the function call operator and both multiplication 

47 and division by a constant via operator overloading. Subclasses provide 

48 those operations by implementing the ``_evaluate`` and 

49 ``_multiply_constant`` hooks (respectively). 

50 

51 This interface will probably change in the future to incorporate options 

52 for dealing with out-of-bounds positions. At present the behavior for 

53 such positions is implementation-specific and should not be relied upon. 

54 """ 

55 

56 @property 

57 @abstractmethod 

58 def bounds(self) -> Bounds: 

59 """The region over which this field can be evaluated (`.Bounds`).""" 

60 raise NotImplementedError() 

61 

62 @property 

63 @abstractmethod 

64 def unit(self) -> astropy.units.UnitBase | None: 

65 """The units of the field (`astropy.units.UnitBase` or `None`).""" 

66 raise NotImplementedError() 

67 

68 @property 

69 @abstractmethod 

70 def is_constant(self) -> bool: 

71 """Whether the field is spatially constant (`bool`).""" 

72 raise NotImplementedError() 

73 

74 @overload 

75 def __call__( 

76 self, point: XY[int | float] | YX[int | float], /, *, quantity: Literal[False] = False 

77 ) -> float: ... 

78 

79 @overload 

80 def __call__( 

81 self, point: XY[int | float] | YX[int | float], /, *, quantity: Literal[True] 

82 ) -> astropy.units.Quantity: ... 

83 

84 @overload 

85 def __call__( 

86 self, point: XY[npt.ArrayLike] | YX[npt.ArrayLike], /, *, quantity: bool = ... 

87 ) -> np.ndarray | astropy.units.Quantity: ... 

88 

89 @overload 

90 def __call__(self, /, *, x: int | float, y: int | float, quantity: Literal[False] = False) -> float: ... 

91 

92 @overload 

93 def __call__( 

94 self, /, *, x: int | float, y: int | float, quantity: Literal[True] 

95 ) -> astropy.units.Quantity: ... 

96 

97 @overload 

98 def __call__( 

99 self, /, *, x: npt.ArrayLike, y: npt.ArrayLike, quantity: bool = ... 

100 ) -> np.ndarray | astropy.units.Quantity: ... 

101 

102 def __call__( 

103 self, 

104 point: XY[Any] | YX[Any] | None = None, 

105 /, 

106 *, 

107 x: Any = None, 

108 y: Any = None, 

109 quantity: bool = False, 

110 ) -> float | np.ndarray | astropy.units.Quantity: 

111 match point: 

112 case None: 

113 if x is None or y is None: 113 ↛ 114line 113 didn't jump to line 114 because the condition on line 113 was never true

114 raise TypeError("Pass either a point or both x= and y= to field call.") 

115 case XY() | YX(): 115 ↛ 119line 115 didn't jump to line 119 because the pattern on line 115 always matched

116 if x is not None or y is not None: 

117 raise TypeError("Field call point argument is mutually exclusive with x= and y=.") 

118 x, y = point.x, point.y 

119 case _: 

120 raise TypeError(f"Unexpected positional argument type: {type(point)!r}.") 

121 x = np.asarray(x) 

122 y = np.asarray(y) 

123 scalar = not np.broadcast(x, y).shape 

124 result = self._evaluate(x=x, y=y, quantity=quantity) 

125 if scalar: 

126 if quantity: 

127 return result # 0-d Quantity 

128 return float(result) 

129 return result 

130 

131 @abstractmethod 

132 def render( 

133 self, 

134 bbox: Box | None = None, 

135 *, 

136 dtype: np.typing.DTypeLike | None = None, 

137 ) -> Image: 

138 """Create an image realization of the field. 

139 

140 Parameters 

141 ---------- 

142 bbox 

143 Bounding box of the image. If not provided, ``self.bounds.bbox`` 

144 will be used. 

145 dtype 

146 Pixel data type for the returned image. 

147 """ 

148 raise NotImplementedError() 

149 

150 def __mul__(self, factor: float | astropy.units.Quantity | astropy.units.UnitBase) -> Self: 

151 return self._multiply_constant(factor) 

152 

153 def __rmul__(self, factor: float | astropy.units.Quantity | astropy.units.UnitBase) -> Self: 

154 return self._multiply_constant(factor) 

155 

156 def __truediv__(self, factor: float | astropy.units.Quantity | astropy.units.UnitBase) -> Self: 

157 return self._multiply_constant(1.0 / factor) 

158 

159 @abstractmethod 

160 def _evaluate( 

161 self, *, x: np.ndarray, y: np.ndarray, quantity: bool 

162 ) -> np.ndarray | astropy.units.Quantity: 

163 """Evaluate at non-gridded points. 

164 

165 Parameters 

166 ---------- 

167 x 

168 X coordinates to evaluate at. 

169 y 

170 Y coordinates to evaluate at; must be broadcast-compatible with 

171 ``x``. 

172 quantity 

173 If `True`, return an `astropy.units.Quantity` instead of a 

174 `numpy.ndarray`. If `unit` is `None`, the returned object will 

175 be a dimensionless `~astropy.units.Quantity`. 

176 """ 

177 raise NotImplementedError() 

178 

179 @abstractmethod 

180 def _multiply_constant(self, factor: float | astropy.units.Quantity | astropy.units.UnitBase) -> Self: 

181 """Multiply by a constant, returning a new field of the same type. 

182 

183 Parameters 

184 ---------- 

185 factor 

186 Factor to multiply by. When this has units, those should multiply 

187 ``self.unit`` or set the units of the returned field if 

188 ``self.unit is None``. 

189 """ 

190 raise NotImplementedError() 

191 

192 def _describe(self, options: DescribeOptions = DescribeOptions(), /) -> Report: 

193 """Return a `Report` describing this field. 

194 

195 Parameters 

196 ---------- 

197 options : `DescribeOptions`, optional 

198 Unused; accepted for interface compatibility. 

199 """ 

200 return Report( 

201 type_name=type(self).__name__, 

202 summary=f"{type(self).__name__} over {self.bounds}", 

203 fields=[ 

204 ReportField(label="bounds", value=self.bounds, role=FieldRole.DERIVED), 

205 ReportField(label="unit", value=self.unit, role=FieldRole.DERIVED), 

206 ReportField(label="is_constant", value=self.is_constant, role=FieldRole.DERIVED), 

207 ], 

208 ) 

209 

210 def to_legacy(self) -> LegacyBoundedField: 

211 """Convert to a legacy `lsst.afw.math.BoundedField`.""" 

212 raise NotImplementedError(f"{type(self).__name__} has no lsst.afw.math.BoundedField representation.") 

213 

214 def to_legacy_photo_calib(self, image_unit: astropy.units.UnitBase) -> LegacyPhotoCalib: 

215 """Convert to a legacy `lsst.afw.image.PhotoCalib`. 

216 

217 Parameters 

218 ---------- 

219 image_unit 

220 The units of the pixels the returned ``PhotoCalib`` will be 

221 associated with. 

222 """ 

223 from lsst.afw.image import PhotoCalib 

224 

225 if (result := self.make_legacy_photo_calib(image_unit)) is not None: 

226 return result 

227 field = self 

228 factor = image_unit.to(astropy.units.nJy / self.unit) 

229 if factor != 1.0: 

230 # TODO[DM-54556]: make sure this shouldn't be 1/factor. 

231 field = self._multiply_constant(factor) # this lies about units, but we'll discard them anyway. 

232 (field_at_center,) = field( 

233 x=np.array([field.bounds.bbox.x.center]), 

234 y=np.array([field.bounds.bbox.y.center]), 

235 ) 

236 if field.is_constant: 

237 return PhotoCalib(field_at_center) 

238 else: 

239 # Constructing a legacy PhotoCalib from a BoundedField alone 

240 # doesn't always work, because ProductBoundedField doesn't 

241 # implement computeMean(). Luckily PhotoCalib doesn't really care 

242 # about getting a true mean; it just wants some sort of central 

243 # tendency, so we can evaluate the field at the bbox center and use 

244 # that (this is what fgcmcal does when it makes a 

245 # ProductBoundedField PhotoCalib). 

246 return PhotoCalib( 

247 calibrationMean=field_at_center, 

248 calibrationErr=0.0, # we don't round-trip this; it's not useful 

249 calibration=field.to_legacy(), 

250 isConstant=False, 

251 ) 

252 

253 @staticmethod 

254 def make_legacy_photo_calib(image_unit: astropy.units.UnitBase) -> LegacyPhotoCalib | None: 

255 """Make a legacy `lsst.afw.image.PhotoCalib` for an image with the 

256 given units, if that is possible without a photometric scaling field. 

257 

258 Parameters 

259 ---------- 

260 image_unit 

261 Units of the image the photometric calibration applies to. 

262 """ 

263 from lsst.afw.image import PhotoCalib 

264 

265 try: 

266 factor = image_unit.to(astropy.units.nJy) 

267 except astropy.units.UnitConversionError: 

268 pass 

269 else: 

270 return PhotoCalib(factor) 

271 return None 

272 

273 def _handle_factor_units( 

274 self, factor: float | astropy.units.Quantity | astropy.units.UnitBase 

275 ) -> tuple[float, astropy.units.UnitBase | None]: 

276 """Interpret the ``factor`` argument to `multiply_constant` and apply 

277 any units it carries to this field's units. 

278 

279 This is a convenience function for subclass implementations of 

280 `multiply_constant`. 

281 

282 Parameters 

283 ---------- 

284 factor 

285 Factor passed by the caller. 

286 

287 Returns 

288 ------- 

289 `float` 

290 The factor to multiply by as a pure `float` 

291 `astropy.units.UnitBase` | `None` 

292 The units for the new field returned by `multiply_constant`. 

293 """ 

294 unit = self.unit 

295 factor_unit = None 

296 if isinstance(factor, astropy.units.Quantity): 

297 factor_unit = factor.unit 

298 factor = factor.to_value() 

299 elif isinstance(factor, astropy.units.UnitBase): 

300 factor_unit = factor 

301 factor = 1.0 

302 if factor_unit is not None: 

303 if unit is None: 

304 unit = factor_unit 

305 else: 

306 unit *= factor_unit 

307 return factor, unit 

308 

309 

310if TYPE_CHECKING: 

311 

312 def _test_types() -> None: 

313 field = cast(BaseField, None) 

314 arr = np.zeros(3) 

315 

316 # Scalar inputs without quantity → float 

317 assert_type(field(x=1.0, y=2.0), float) 

318 assert_type(field(x=1.0, y=2.0, quantity=False), float) 

319 

320 # Scalar inputs with quantity=True → astropy.units.Quantity 

321 assert_type(field(x=1.0, y=2.0, quantity=True), astropy.units.Quantity) 

322 

323 # Array-like inputs → np.ndarray | astropy.units.Quantity 

324 assert_type(field(x=arr, y=arr), np.ndarray | astropy.units.Quantity)