Skip to content

Commit

Permalink
Merge pull request #73 from lsst/tickets/DM-47043
Browse files Browse the repository at this point in the history
DM-47043: Fix bug with healpix envelope related to pixels with incomplete neighbors
  • Loading branch information
erykoff authored Oct 29, 2024
2 parents fc885e0 + e18d690 commit 781aacb
Show file tree
Hide file tree
Showing 2 changed files with 28 additions and 2 deletions.
4 changes: 2 additions & 2 deletions python/lsst/sphgeom/_healpixPixelization.py
Original file line number Diff line number Diff line change
Expand Up @@ -89,9 +89,9 @@ def envelope(self, region: Region, maxRanges: int = 0):

# Dilate the high resolution pixels by one to ensure that the full
# region is completely covered at high resolution.
neighbors = hpg.neighbors(self._nside_highres, pixels_highres)
neighbors = hpg.neighbors(self._nside_highres, pixels_highres).ravel()
# Shift back to the original resolution and uniquify
pixels = np.unique(np.right_shift(neighbors.ravel(), self._bit_shift))
pixels = np.unique(np.right_shift(neighbors[neighbors >= 0], self._bit_shift))

return RangeSet(pixels)

Expand Down
26 changes: 26 additions & 0 deletions tests/test_HealpixPixelization.py
Original file line number Diff line number Diff line change
Expand Up @@ -114,6 +114,32 @@ def test_envelope(self):
)
self._check_envelope(h, circle, [98, 99, 104, 105])

def test_envelope_missing_neighbors(self):
"""Test envelope, with a pixel with missing neighbors.
Testing DM-47043.
"""
h = HealpixPixelization(0)
h_highres = HealpixPixelization(4)

self.assertEqual(h._nside_highres, h_highres._nside)

# Find a pixel with a missing neighbor at high res.
n = hpg.neighbors(h_highres._nside, np.arange(hpg.nside_to_npixel(h_highres._nside)))
ind = np.argmin(n)
self.assertEqual(n.ravel()[ind], -1)

# This is the pixel with a bad (-1) neighbor.
badpix = ind // 8
ra, dec = hpg.pixel_to_angle(h_highres._nside, badpix)

box = Box(
point1=LonLat.fromDegrees(ra - 0.1, dec - 0.1),
point2=LonLat.fromDegrees(ra + 0.1, dec + 0.1),
)

self._check_envelope(h, box, [0, 1, 5])

def _check_envelope(self, pixelization, region, check_pixels):
"""Check the envelope from a region.
Expand Down

0 comments on commit 781aacb

Please sign in to comment.