Coverage for python / lsst / images / _cell_grid.py: 44%

109 statements  

« prev     ^ index     » next       coverage.py v7.14.0, created at 2026-05-16 07:54 +0000

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# This module is conceptually part of the 'cells' subpackage, but we don't 

15# want the stuff in '_concrete_bounds' to depend on all of that. So the 

16# basic CellGrid and CellGridBounds objects are defined here, used in both 

17# places, and exported from 'cells'. 

18 

19__all__ = ( 

20 "CellGrid", 

21 "CellGridBounds", 

22 "CellIJ", 

23 "PatchDefinition", 

24) 

25 

26import dataclasses 

27import math 

28from collections.abc import Iterator 

29from functools import cached_property 

30from typing import TYPE_CHECKING, Any, overload 

31 

32import numpy as np 

33import pydantic 

34 

35from ._geom import YX, Bounds, Box 

36 

37if TYPE_CHECKING: 

38 try: 

39 from lsst.cell_coadds import UniformGrid 

40 from lsst.skymap import Index2D 

41 except ImportError: 

42 type UniformGrid = Any # type: ignore[no-redef] 

43 type Index2D = Any # type: ignore[no-redef] 

44 

45 

46@dataclasses.dataclass(frozen=True, order=True) 

47class CellIJ: 

48 """An index in a grid of cells. 

49 

50 Notes 

51 ----- 

52 This is deliberately not a `tuple` or other `~collections.abc.Sequence` in 

53 order to make it typing-incompatible with sequence-based pixel coordinate 

54 pairs (e.g. `.YX`). This also allows it to have addition and subtraction 

55 operators. 

56 """ 

57 

58 i: int 

59 """The y / row object.""" 

60 

61 j: int 

62 """The x / column object.""" 

63 

64 def __add__(self, other: CellIJ) -> CellIJ: 

65 return CellIJ(i=self.i + other.i, j=self.j + other.j) 

66 

67 def __sub__(self, other: CellIJ) -> CellIJ: 

68 return CellIJ(i=self.i - other.i, j=self.j - other.j) 

69 

70 @staticmethod 

71 def from_legacy(legacy_index: Index2D) -> CellIJ: 

72 """Convert from a legacy `lsst.skymap.Index2D` instance. 

73 

74 Notes 

75 ----- 

76 `lsst.skymap.Index2D` is ordered ``(x, y)``, i.e. ``(j, i)``. 

77 """ 

78 return CellIJ(i=legacy_index.y, j=legacy_index.x) 

79 

80 def to_legacy(self) -> Index2D: 

81 """Convert to a legacy `lsst.skymap.Index2D` instance. 

82 

83 Notes 

84 ----- 

85 `lsst.skymap.Index2D` is ordered ``(x, y)``, i.e. ``(j, i)``. 

86 """ 

87 from lsst.skymap import Index2D 

88 

89 return Index2D(x=self.j, y=self.i) 

90 

91 

92class CellGrid(pydantic.BaseModel, frozen=True): 

93 """A grid of rectangular cells with no overlaps or space between cells. 

94 

95 Notes 

96 ----- 

97 A cell grid usually corresponds to a full patch, but we do not explicitly 

98 encode this in the type to permit full-tract grids, which would have to 

99 drop the cells in patch overlap regions and re-label all cells. 

100 

101 Subsets of grids are usually represented via `CellGridBounds`. 

102 """ 

103 

104 bbox: Box = pydantic.Field( 

105 description=( 

106 "Bounding box of the grid of cells (snapped to cell boundaries. " 

107 "The cell with index (i=0, j=0) always has a corner at ``(y=bbox.y.min, x=bbox.x.min)`` " 

108 "but there is no expectation that ``(y=bbox.y.min, x=bbox.x.min)`` be ``(y=0, x=0)``." 

109 ) 

110 ) 

111 cell_shape: YX[int] = pydantic.Field(description="Shape of each cell in pixels.") 

112 

113 @property 

114 def grid_shape(self) -> CellIJ: 

115 """The number of cells in each dimension (`CellIJ`).""" 

116 return CellIJ(i=self.bbox.y.size // self.cell_shape.y, j=self.bbox.x.size // self.cell_shape.x) 

117 

118 def index_of(self, *, y: int, x: int) -> CellIJ: 

119 """Return the 2-d index of the cell that contains the given pixel. 

120 

121 Parameters 

122 ---------- 

123 y 

124 Y cell index. 

125 x 

126 X cell index. 

127 """ 

128 return CellIJ( 

129 i=(y - self.bbox.y.start) // self.cell_shape.y, 

130 j=(x - self.bbox.x.start) // self.cell_shape.x, 

131 ) 

132 

133 def bbox_of(self, cell: CellIJ) -> Box: 

134 """Return the bounding box of the given cell.""" 

135 return Box.from_shape( 

136 self.cell_shape, 

137 start=YX( 

138 y=cell.i * self.cell_shape.y + self.bbox.y.start, 

139 x=cell.j * self.cell_shape.x + self.bbox.x.start, 

140 ), 

141 ) 

142 

143 @staticmethod 

144 def from_legacy(legacy: UniformGrid) -> CellGrid: 

145 """Construct from a legacy grid object. 

146 

147 Parameters 

148 ---------- 

149 legacy 

150 Legacy grid to convert. 

151 """ 

152 if legacy.padding: 

153 raise ValueError("Only cell grids with no padding are supported.") 

154 bbox = Box.from_legacy(legacy.bbox) 

155 cell_shape = YX(y=legacy.cell_size.y, x=legacy.cell_size.x) 

156 return CellGrid(bbox=bbox, cell_shape=cell_shape) 

157 

158 

159class CellGridBounds(pydantic.BaseModel, frozen=True): 

160 """A region of pixels defined by a set of cells within a grid. 

161 

162 Notes 

163 ----- 

164 This data structure is optimized for the case where a continguous 

165 rectangular region of the grid (the `bbox` attribute) is populated with 

166 only a few exceptions (the `missing` set). 

167 

168 Slicing a `CellGridBounds` with a `.Box` returns a new `CellGridBounds` 

169 with just the cells that overlap that box. As always, 

170 `CellGridBounds.bbox` will be snapped to the outer boundaries of those 

171 cells, so it will contain (and not generally equal) the given box. 

172 """ 

173 

174 grid: CellGrid = pydantic.Field(description="Definition of the grid that defines the cells.") 

175 bbox: Box = pydantic.Field(description="Pixel bounding box of the region (snapped to cell boundaries).") 

176 missing: frozenset[CellIJ] = pydantic.Field( 

177 default=frozenset(), 

178 description=( 

179 "Indices of cells that are missing, where (i=0, j=0) is the cell that starts at grid.bbox.start." 

180 ), 

181 ) 

182 

183 @cached_property 

184 def grid_start(self) -> CellIJ: 

185 """The index of the first cell in this bounds' bounding box within 

186 its grid. 

187 """ 

188 return self.grid.index_of(y=self.bbox.y.start, x=self.bbox.x.start) 

189 

190 @cached_property 

191 def grid_stop(self) -> CellIJ: 

192 """One-past-the-last indices for the cells in these bounds, within 

193 its grid. 

194 """ 

195 return self.grid.index_of(y=self.bbox.y.stop, x=self.bbox.x.stop) 

196 

197 @overload 

198 def contains(self, *, x: int, y: int) -> bool: ... 198 ↛ exitline 198 didn't return from function 'contains' because

199 

200 @overload 

201 def contains(self, *, x: np.ndarray, y: np.ndarray) -> np.ndarray: ... 201 ↛ exitline 201 didn't return from function 'contains' because

202 

203 def contains(self, *, x: Any, y: Any) -> Any: 

204 """Test whether these bounds contain one or more points. 

205 

206 Parameters 

207 ---------- 

208 x 

209 One or more integer X coordinates to test for containment. 

210 If an array, an array of results will be returned. 

211 y 

212 One or more integer Y coordinates to test for containment. 

213 If an array, an array of results will be returned. 

214 

215 Returns 

216 ------- 

217 `bool` | `numpy.ndarray` 

218 If ``x`` and ``y`` are both scalars, a single `bool` value. If 

219 ``x`` and ``y`` are arrays, a boolean array with their broadcasted 

220 shape. 

221 """ 

222 result = self.bbox.contains(x=x, y=y) 

223 if not self.missing: 

224 return result 

225 match result: 

226 case False: 

227 return False 

228 case True: 

229 return self.grid.index_of(x=x, y=y) not in self.missing 

230 case np.ndarray(): 

231 for box in self.missing_boxes(): 

232 result = np.logical_and(result, np.logical_not(box.contains(x=x, y=y))) 

233 return result 

234 

235 def intersection(self, other: Bounds) -> Bounds: 

236 """Compute the intersection of this bounds object with another.""" 

237 from ._concrete_bounds import _intersect_cgb 

238 

239 return _intersect_cgb(self, other) 

240 

241 def contains_cell(self, index: CellIJ) -> bool: 

242 """Test whether the given cell is in the bounds.""" 

243 return ( 

244 (index.i >= self.grid_start.i and index.i < self.grid_stop.i) 

245 and (index.j >= self.grid_start.j and index.j < self.grid_stop.j) 

246 and index not in self.missing 

247 ) 

248 

249 def missing_boxes(self) -> Iterator[Box]: 

250 """Iterate over the bounding boxes of the missing cells.""" 

251 for index in sorted(self.missing): 

252 yield self.grid.bbox_of(index) 

253 

254 def cell_indices(self) -> Iterator[CellIJ]: 

255 """Iterate over the indices of the cells in these bounds.""" 

256 for i in range(self.grid_start.i, self.grid_stop.i): 

257 for j in range(self.grid_start.j, self.grid_stop.j): 

258 index = CellIJ(i=i, j=j) 

259 if index not in self.missing: 

260 yield index 

261 

262 def __getitem__(self, bbox: Box) -> CellGridBounds: 

263 if not self.bbox.contains(bbox): 

264 raise ValueError( 

265 f"Original grid bounding box {self.bbox} does not contain the subset bounding box {bbox}." 

266 ) 

267 c = self.grid.cell_shape 

268 s = self.grid.bbox.start 

269 i1 = (bbox.y.start - s.y) // c.y 

270 j1 = (bbox.x.start - s.x) // c.x 

271 i2 = math.ceil((bbox.y.stop - s.y) / c.y) 

272 j2 = math.ceil((bbox.x.stop - s.x) / c.x) 

273 subset_bbox = Box.factory[i1 * c.y + s.y : i2 * c.y + s.y, j1 * c.x + s.x : j2 * c.x + s.x] 

274 grid_as_box = Box.factory[i1:i2, j1:j2] 

275 subset_missing = {index for index in self.missing if grid_as_box.contains(y=index.i, x=index.j)} 

276 return CellGridBounds(grid=self.grid, bbox=subset_bbox, missing=frozenset(subset_missing)) 

277 

278 def serialize(self) -> CellGridBounds: 

279 """Convert a bounds instance into a serializable object.""" 

280 return self 

281 

282 def deserialize(self) -> CellGridBounds: 

283 """Deserialize a bounds object on the assumption it is a 

284 `CellGridBounds`. 

285 

286 This method just returns the `CellGridBounds` itself, since that 

287 already provides Pydantic serialization hooks. It exists for 

288 compatibility with the `.Bounds` protocol. 

289 """ 

290 return self 

291 

292 

293class PatchDefinition(pydantic.BaseModel, frozen=True): 

294 """Identifiers and geometry for a full patch.""" 

295 

296 id: int = pydantic.Field(description="ID for the patch.") 

297 index: YX[int] = pydantic.Field(description="2-d index of this patch within the tract.") 

298 inner_bbox: Box = pydantic.Field(description="Inner bounding box of this patch.") 

299 cells: CellGrid = pydantic.Field(description="Cell grid for the full patch.") 

300 

301 @property 

302 def outer_bbox(self) -> Box: 

303 """The outer bounding box of this patch (`.Box`).""" 

304 return self.cells.bbox