Coverage for tests/test_stcs.py: 17%
128 statements
« prev ^ index » next coverage.py v7.14.1, created at 2026-05-29 01:11 -0700
1# This file is part of sphgeom.
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 software is dual licensed under the GNU General Public License and also
10# under a 3-clause BSD license. Recipients may choose which of these licenses
11# to use; please see the files gpl-3.0.txt and/or bsd_license.txt,
12# respectively. If you choose the GPL option then the following text applies
13# (but note that there is still no warranty even if you opt for BSD instead):
14#
15# This program is free software: you can redistribute it and/or modify
16# it under the terms of the GNU General Public License as published by
17# the Free Software Foundation, either version 3 of the License, or
18# (at your option) any later version.
19#
20# This program is distributed in the hope that it will be useful,
21# but WITHOUT ANY WARRANTY; without even the implied warranty of
22# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
23# GNU General Public License for more details.
24#
25# You should have received a copy of the GNU General Public License
26# along with this program. If not, see <http://www.gnu.org/licenses/>.
28import math
29import unittest
31from lsst.sphgeom import (
32 Angle,
33 Box,
34 Circle,
35 ConvexPolygon,
36 Ellipse,
37 IntersectionRegion,
38 LonLat,
39 Region,
40 UnionRegion,
41 UnitVector3d,
42)
45class StcsTestCase(unittest.TestCase):
46 """Test STC-S string generation."""
48 def assert_stcs_equal(self, stcs1: str, stcs2: str):
49 """Compare two STC-S strings, treating numeric tokens with a
50 floating-point tolerance and string tokens (shape/operator/frame
51 keywords, parentheses) by exact match.
52 """
53 toks1 = stcs1.replace("(", " ( ").replace(")", " ) ").split()
54 toks2 = stcs2.replace("(", " ( ").replace(")", " ) ").split()
55 self.assertEqual(len(toks1), len(toks2), f"{stcs1!r} vs {stcs2!r}")
56 for t1, t2 in zip(toks1, toks2, strict=True):
57 try:
58 f1 = float(t1)
59 f2 = float(t2)
60 except ValueError:
61 self.assertEqual(t1, t2)
62 else:
63 self.assertAlmostEqual(f1, f2, places=6)
65 def test_base_region_default_raises(self):
66 """Calling ``to_ivoa_stcs`` through ``Region`` on a subclass that
67 cannot be represented in STC-S raises ``NotImplementedError``.
68 """
69 box = Box(LonLat.fromDegrees(1.0, 2.0), LonLat.fromDegrees(5.0, 6.0))
70 with self.assertRaises(NotImplementedError):
71 Region.to_ivoa_stcs(box)
73 def test_circle(self):
74 """Circle round-trips through STC-S."""
75 circle = Circle(
76 UnitVector3d(LonLat.fromDegrees(180.0, 30.0)),
77 Angle.fromDegrees(2.0),
78 )
79 self.assert_stcs_equal(
80 circle.to_ivoa_stcs(),
81 "Circle ICRS 180.0 30.0 2.0",
82 )
84 def test_circle_frame_argument(self):
85 """The frame argument is emitted verbatim."""
86 circle = Circle(
87 UnitVector3d(LonLat.fromDegrees(180.0, 30.0)),
88 Angle.fromDegrees(2.0),
89 )
90 self.assert_stcs_equal(
91 circle.to_ivoa_stcs(frame="GALACTIC"),
92 "Circle GALACTIC 180.0 30.0 2.0",
93 )
95 def test_circle_no_frame(self):
96 """Passing an empty frame omits the frame keyword."""
97 circle = Circle(
98 UnitVector3d(LonLat.fromDegrees(180.0, 30.0)),
99 Angle.fromDegrees(2.0),
100 )
101 self.assert_stcs_equal(
102 circle.to_ivoa_stcs(frame=""),
103 "Circle 180.0 30.0 2.0",
104 )
106 def test_polygon(self):
107 """Polygon emits Polygon <frame> followed by lon/lat pairs."""
108 vertices = [
109 UnitVector3d(LonLat.fromDegrees(12.0, 34.0)),
110 UnitVector3d(LonLat.fromDegrees(14.0, 34.0)),
111 UnitVector3d(LonLat.fromDegrees(14.0, 36.0)),
112 UnitVector3d(LonLat.fromDegrees(12.0, 36.0)),
113 ]
114 poly = ConvexPolygon(vertices)
115 # ConvexPolygon may rotate the vertex order; compare token counts and
116 # numeric values pairwise after extracting (lon, lat) pairs.
117 stcs = poly.to_ivoa_stcs()
118 tokens = stcs.split()
119 self.assertEqual(tokens[0], "Polygon")
120 self.assertEqual(tokens[1], "ICRS")
121 coords = [float(t) for t in tokens[2:]]
122 self.assertEqual(len(coords), 8)
123 # Recover the 4 (lon, lat) pairs as a set of rounded tuples.
124 emitted = {(round(coords[i], 6), round(coords[i + 1], 6)) for i in range(0, 8, 2)}
125 expected = {(12.0, 34.0), (14.0, 34.0), (14.0, 36.0), (12.0, 36.0)}
126 self.assertEqual(emitted, expected)
128 def test_polygon_roundtrip_precision(self):
129 """Lon/lat values round-trip through STC-S without precision loss.
131 This locks in shortest-round-trip formatting: when the body string is
132 re-parsed with ``float()``, every value must equal the polygon's
133 actual stored vertex coordinate exactly (no tolerance).
134 """
135 vertices = [
136 UnitVector3d(LonLat.fromDegrees(12.345678901234567, 34.56789012345678)),
137 UnitVector3d(LonLat.fromDegrees(56.789012345678901, 78.90123456789012)),
138 UnitVector3d(LonLat.fromDegrees(-23.456789012345678, -45.67890123456789)),
139 ]
140 poly = ConvexPolygon(vertices)
141 stcs = poly.to_ivoa_stcs()
142 tokens = stcs.split()
143 self.assertEqual(tokens[0], "Polygon")
144 self.assertEqual(tokens[1], "ICRS")
145 emitted = [float(t) for t in tokens[2:]]
146 # Compare to the polygon's stored vertices (ConvexPolygon may reorder
147 # or rotate the input, so compare against its actual getVertices()).
148 expected = []
149 for v in poly.getVertices():
150 ll = LonLat(v)
151 expected.append(ll.getLon().asDegrees())
152 expected.append(ll.getLat().asDegrees())
153 self.assertEqual(len(emitted), len(expected))
154 for got, want in zip(emitted, expected, strict=True):
155 self.assertEqual(got, want, f"got {got!r} expected {want!r}")
157 def test_ellipse(self):
158 """Ellipse round-trips through STC-S, including position angle."""
159 center = UnitVector3d(LonLat.fromDegrees(180.0, 30.0))
160 ellipse = Ellipse(
161 center,
162 Angle.fromDegrees(2.0), # alpha
163 Angle.fromDegrees(1.0), # beta
164 Angle.fromDegrees(45.0), # orientation (PA)
165 )
166 self.assert_stcs_equal(
167 ellipse.to_ivoa_stcs(),
168 "Ellipse ICRS 180.0 30.0 2.0 1.0 45.0",
169 )
171 def test_ellipse_pa_zero(self):
172 """Ellipse with orientation=0 emits PA close to 0."""
173 center = UnitVector3d(LonLat.fromDegrees(45.0, 0.0))
174 ellipse = Ellipse(
175 center,
176 Angle.fromDegrees(2.0),
177 Angle.fromDegrees(1.0),
178 Angle.fromDegrees(0.0),
179 )
180 tokens = ellipse.to_ivoa_stcs().split()
181 # Last token is the position angle.
182 pa = float(tokens[-1])
183 # PA may come back as 0 or 180 (axis is unsigned).
184 self.assertTrue(
185 math.isclose(pa, 0.0, abs_tol=1e-6) or math.isclose(pa, 180.0, abs_tol=1e-6),
186 f"unexpected PA {pa}",
187 )
189 def test_ellipse_no_frame(self):
190 """Passing an empty frame omits the frame keyword."""
191 center = UnitVector3d(LonLat.fromDegrees(180.0, 30.0))
192 ellipse = Ellipse(
193 center,
194 Angle.fromDegrees(2.0),
195 Angle.fromDegrees(1.0),
196 Angle.fromDegrees(45.0),
197 )
198 self.assert_stcs_equal(
199 ellipse.to_ivoa_stcs(frame=""),
200 "Ellipse 180.0 30.0 2.0 1.0 45.0",
201 )
203 def test_ellipse_full_raises(self):
204 """A full ellipse cannot be represented as STC-S."""
205 center = UnitVector3d(LonLat.fromDegrees(0.0, 0.0))
206 ellipse = Ellipse(center, Angle.fromDegrees(180.0)) # full
207 self.assertTrue(ellipse.isFull())
208 with self.assertRaises(ValueError):
209 ellipse.to_ivoa_stcs()
211 def test_ellipse_empty_raises(self):
212 """An empty ellipse cannot be represented as STC-S."""
213 # Default-constructed Ellipse is empty.
214 ellipse = Ellipse()
215 self.assertTrue(ellipse.isEmpty())
216 with self.assertRaises(ValueError):
217 ellipse.to_ivoa_stcs()
219 def test_box_not_supported(self):
220 """Box explicitly raises NotImplementedError with a helpful message."""
221 box = Box(
222 LonLat.fromDegrees(1.0, 2.0),
223 LonLat.fromDegrees(5.0, 6.0),
224 )
225 with self.assertRaises(NotImplementedError) as cm:
226 box.to_ivoa_stcs()
227 # The error message should mention an alternative for callers.
228 self.assertIn("Polygon", str(cm.exception))
230 def test_union(self):
231 """Union of two circles emits a single ICRS keyword."""
232 c1 = Circle(UnitVector3d(LonLat.fromDegrees(180.0, 10.0)), Angle.fromDegrees(2.0))
233 c2 = Circle(UnitVector3d(LonLat.fromDegrees(190.0, 20.0)), Angle.fromDegrees(1.0))
234 u = UnionRegion(c1, c2)
235 self.assert_stcs_equal(
236 u.to_ivoa_stcs(),
237 "Union ICRS ( Circle 180.0 10.0 2.0 Circle 190.0 20.0 1.0 )",
238 )
240 def test_union_nested(self):
241 """Nested unions are flattened by ``UnionRegion`` at construction
242 time, so the resulting STC-S is a single ``Union`` with all three
243 operands and a single frame keyword.
244 """
245 c1 = Circle(UnitVector3d(LonLat.fromDegrees(0.0, 0.0)), Angle.fromDegrees(1.0))
246 c2 = Circle(UnitVector3d(LonLat.fromDegrees(10.0, 0.0)), Angle.fromDegrees(1.0))
247 c3 = Circle(UnitVector3d(LonLat.fromDegrees(20.0, 0.0)), Angle.fromDegrees(1.0))
248 nested = UnionRegion(c1, UnionRegion(c2, c3))
249 # UnionRegion flattens same-kind operands automatically.
250 self.assertEqual(nested.nOperands(), 3)
251 stcs = nested.to_ivoa_stcs()
252 # ICRS appears exactly once.
253 self.assertEqual(stcs.count("ICRS"), 1)
254 # Only one Union keyword (flattened).
255 self.assertEqual(stcs.count("Union"), 1)
256 self.assert_stcs_equal(
257 stcs,
258 "Union ICRS ( Circle 0.0 0.0 1.0 Circle 10.0 0.0 1.0 Circle 20.0 0.0 1.0 )",
259 )
261 def test_union_no_frame(self):
262 """Passing an empty frame omits the frame keyword at the top level."""
263 c1 = Circle(UnitVector3d(LonLat.fromDegrees(180.0, 10.0)), Angle.fromDegrees(2.0))
264 c2 = Circle(UnitVector3d(LonLat.fromDegrees(190.0, 20.0)), Angle.fromDegrees(1.0))
265 u = UnionRegion(c1, c2)
266 self.assert_stcs_equal(
267 u.to_ivoa_stcs(frame=""),
268 "Union ( Circle 180.0 10.0 2.0 Circle 190.0 20.0 1.0 )",
269 )
271 def test_union_with_unsupported_operand(self):
272 """A Union containing a Box raises NotImplementedError."""
273 circle = Circle(UnitVector3d(LonLat.fromDegrees(0.0, 0.0)), Angle.fromDegrees(1.0))
274 box = Box(LonLat.fromDegrees(1.0, 2.0), LonLat.fromDegrees(5.0, 6.0))
275 u = UnionRegion(circle, box)
276 with self.assertRaises(NotImplementedError):
277 u.to_ivoa_stcs()
279 def test_intersection(self):
280 """Intersection of a circle and a polygon."""
281 circle = Circle(UnitVector3d(LonLat.fromDegrees(0.0, 0.0)), Angle.fromDegrees(5.0))
282 vertices = [
283 UnitVector3d(LonLat.fromDegrees(-1.0, -1.0)),
284 UnitVector3d(LonLat.fromDegrees(1.0, -1.0)),
285 UnitVector3d(LonLat.fromDegrees(1.0, 1.0)),
286 UnitVector3d(LonLat.fromDegrees(-1.0, 1.0)),
287 ]
288 poly = ConvexPolygon(vertices)
289 inter = IntersectionRegion(circle, poly)
290 stcs = inter.to_ivoa_stcs()
291 # Polygon vertex order is not deterministic (ConvexPolygon may
292 # rotate), so check the Circle portion exactly and only assert
293 # Polygon presence. ``assert_stcs_equal`` handles the numeric
294 # tolerance (e.g. ``5`` from ``std::to_chars`` matching ``5.0``).
295 polygon_idx = stcs.index(" Polygon ")
296 self.assert_stcs_equal(stcs[:polygon_idx], "Intersection ICRS ( Circle 0.0 0.0 5.0")
297 self.assertEqual(stcs.count("ICRS"), 1)
298 self.assertTrue(stcs.endswith(" )"))
300 def test_intersection_with_unsupported_operand(self):
301 """Intersection containing a Box raises NotImplementedError."""
302 circle = Circle(UnitVector3d(LonLat.fromDegrees(0.0, 0.0)), Angle.fromDegrees(1.0))
303 box = Box(LonLat.fromDegrees(1.0, 2.0), LonLat.fromDegrees(5.0, 6.0))
304 inter = IntersectionRegion(circle, box)
305 with self.assertRaises(NotImplementedError):
306 inter.to_ivoa_stcs()
309if __name__ == "__main__":
310 unittest.main()