Coverage for python/lsst/dax/apdb/pixelization.py: 71%

61 statements  

« prev     ^ index     » next       coverage.py v7.16.0, created at 2026-09-26 02:03 -0700

1# This file is part of dax_apdb. 

2# 

3# Developed for the LSST Data Management System. 

4# This product includes software developed by the LSST Project 

5# (http://www.lsst.org). 

6# See the COPYRIGHT file at the top-level directory of this distribution 

7# for details of code ownership. 

8# 

9# This program is free software: you can redistribute it and/or modify 

10# it under the terms of the GNU General Public License as published by 

11# the Free Software Foundation, either version 3 of the License, or 

12# (at your option) any later version. 

13# 

14# This program is distributed in the hope that it will be useful, 

15# but WITHOUT ANY WARRANTY; without even the implied warranty of 

16# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 

17# GNU General Public License for more details. 

18# 

19# You should have received a copy of the GNU General Public License 

20# along with this program. If not, see <http://www.gnu.org/licenses/>. 

21 

22from __future__ import annotations 

23 

24__all__ = ["Pixelization"] 

25 

26import logging 

27from typing import Any, overload 

28 

29import lsst.sphgeom 

30 

31_LOG = logging.getLogger(__name__) 

32 

33 

34class Pixelization: 

35 """Wrapper for pixelization classes from `sphgeom` with configurable 

36 pixelization type and parameters. 

37 

38 Parameters 

39 ---------- 

40 pixelization : `str` 

41 Name of a pixelization type, one of ""htm", "q3c", "mq3c", or 

42 "healpix". 

43 pix_level : `int` 

44 Pixelization level. 

45 pix_max_ranges : `int` 

46 Maximum number of ranges returned from `envelope()` method. 

47 """ 

48 

49 def __init__(self, pixelization: str, pix_level: int, pix_max_ranges: int): 

50 self._pix_max_ranges = pix_max_ranges 

51 self._is_healpix = False 

52 

53 self.pixelator: lsst.sphgeom.Pixelization 

54 self.level = pix_level 

55 if pixelization == "htm": 

56 self.pixelator = lsst.sphgeom.HtmPixelization(pix_level) 

57 elif pixelization == "q3c": 57 ↛ 58line 57 didn't jump to line 58 because the condition on line 57 was never true

58 self.pixelator = lsst.sphgeom.Q3cPixelization(pix_level) 

59 elif pixelization == "mq3c": 59 ↛ 61line 59 didn't jump to line 61 because the condition on line 59 was always true

60 self.pixelator = lsst.sphgeom.Mq3cPixelization(pix_level) 

61 elif pixelization == "healpix": 

62 # Healpix does not support maxRanges. 

63 self._pix_max_ranges = 0 

64 self._is_healpix = True 

65 self.pixelator = lsst.sphgeom.HealpixPixelization(pix_level) # type: ignore[attr-defined] 

66 else: 

67 raise ValueError(f"unknown pixelization: {pixelization}") 

68 

69 def pixels(self, region: lsst.sphgeom.Region) -> list[int]: 

70 """Compute set of the pixel indices for given region. 

71 

72 Parameters 

73 ---------- 

74 region : `lsst.sphgeom.Region` 

75 """ 

76 # We want finest set of pixels, so ask as many pixel as reasonable, but 

77 # healpix does not support non-zero maxRanges. 

78 ranges = self.pixelator.envelope(region, 0 if self._is_healpix else 1_000_000) 

79 indices = [] 

80 for lower, upper in ranges: 

81 indices += list(range(lower, upper)) 

82 return indices 

83 

84 def circle_pixels(self, ra: float, dec: float, pad_arcsec: float) -> list[int]: 

85 """Make a list of spatial partitions that a small circle touches. 

86 

87 Parameters 

88 ---------- 

89 ra, dec : `float` 

90 Center of a circle, degrees. 

91 pad_arcsec : `float` 

92 Radius of a circle in arcseconds. 

93 

94 Returns 

95 ------- 

96 pixels : `list` [`int`] 

97 All pixels that envelop the circle. 

98 """ 

99 lon_lat = lsst.sphgeom.LonLat.fromDegrees(ra, dec) 

100 center = lsst.sphgeom.UnitVector3d(lon_lat) 

101 region = lsst.sphgeom.Circle(center, lsst.sphgeom.Angle.fromDegrees(pad_arcsec / 3600.0)) 

102 return self.pixels(region) 

103 

104 @overload 

105 def pixel(self, direction: lsst.sphgeom.UnitVector3d, /) -> int: ... 105 ↛ exitline 105 didn't return from function 'pixel' because

106 

107 @overload 

108 def pixel(self, ra: float, dec: float, /) -> int: ... 108 ↛ exitline 108 didn't return from function 'pixel' because

109 

110 def pixel(self, *args: Any) -> int: 

111 """Compute the index of the pixel for given direction. 

112 

113 Parameters 

114 ---------- 

115 args 

116 The method can take either a single `lsst.sphgeom.UnitVector3d` or 

117 a pair of floating point numbers (or values convertible to floats) 

118 representing RA and Dec in degrees. 

119 

120 Returns 

121 ------- 

122 pixel : `int` 

123 Pixel index. 

124 """ 

125 match args: 

126 case (lsst.sphgeom.UnitVector3d() as direction,): 

127 pass 

128 case (ra, dec): 128 ↛ 135line 128 didn't jump to line 135 because the pattern on line 128 always matched

129 try: 

130 direction = lsst.sphgeom.UnitVector3d( 

131 lsst.sphgeom.LonLat.fromDegrees(float(ra), float(dec)) 

132 ) 

133 except (TypeError, ValueError) as exc: 

134 raise TypeError(f"Unexpected arguments: {args}") from exc 

135 case _: 

136 raise TypeError(f"Unexpected arguments: {args}") 

137 index = self.pixelator.index(direction) 

138 return index 

139 

140 def region(self, pixel: int) -> lsst.sphgeom.Region: 

141 """Return region corresponding to a pixel index. 

142 

143 Parameters 

144 ---------- 

145 pixel : `int` 

146 Pixel index. 

147 

148 Returns 

149 ------- 

150 region : `lsst.sphgeom.Region` 

151 Region for a given pixel index. 

152 """ 

153 region = self.pixelator.pixel(pixel) 

154 return region 

155 

156 def envelope(self, region: lsst.sphgeom.Region) -> list[tuple[int, int]]: 

157 """Generate a set of HTM indices covering specified region. 

158 

159 Parameters 

160 ---------- 

161 region: `sphgeom.Region` 

162 Region that needs to be indexed. 

163 

164 Returns 

165 ------- 

166 ranges : `list` of `tuple` 

167 Sequence of ranges, range is a tuple (minHtmID, maxHtmID). 

168 """ 

169 _LOG.debug("region: %s", region) 

170 indices = self.pixelator.envelope(region, self._pix_max_ranges) 

171 

172 if _LOG.isEnabledFor(logging.DEBUG): 172 ↛ 180line 172 didn't jump to line 180 because the condition on line 172 was always true

173 for irange in indices.ranges(): 

174 _LOG.debug( 

175 "range: %s %s", 

176 self.pixelator.toString(irange[0]), 

177 self.pixelator.toString(irange[1]), 

178 ) 

179 

180 return indices.ranges()