fix healpy mask growth to work beyond small angle approx

lsst · Oct 8, 2024 · 88702db · 88702db
1 parent 5d27a22
commit 88702db
Show file tree

Hide file tree

Showing 2 changed files with 32 additions and 16 deletions.
diff --git a/rubin_scheduler/utils/healpy_utils.py b/rubin_scheduler/utils/healpy_utils.py
@@ -398,14 +398,20 @@ def _hp_grow_mask(nside, masked_indx_tuple, grow_size=np.radians(2.0), scale=100
     # get the points that way and save an _hpid2_ra_dec call.
     ra, dec = _hpid2_ra_dec(nside, np.arange(hp.nside2npix(nside))[list(masked_indx_tuple)])
     x, y, z = _xyz_from_ra_dec(ra, dec)
+
+    # convert angular distance to Euclidian distance
+    zp = _xyz_from_ra_dec(0, 0)
+    new_pt = _xyz_from_ra_dec(0, grow_size)
+    dist = np.linalg.norm(zp - new_pt)
+
     if scale is not None:
         x = np.round(x * scale).astype(int)
         y = np.round(y * scale).astype(int)
         z = np.round(z * scale).astype(int)
-        grow_size = np.round(grow_size * scale).astype(int)
+        dist = np.round(dist * scale).astype(int)
     # If there are lots of points, may want to use query_ball_tree
     # instead for speed.
-    lists_of_neighbors = _hp_grow_mask.tree.query_ball_point(np.vstack([x, y, z]).T, grow_size)
+    lists_of_neighbors = _hp_grow_mask.tree.query_ball_point(np.vstack([x, y, z]).T, dist)
 
     u_indx = np.unique(np.concatenate(lists_of_neighbors))
 

diff --git a/tests/utils/test_healutils.py b/tests/utils/test_healutils.py
@@ -102,24 +102,34 @@ def test_mask_grow(self):
         # Need to turn off the machine precision kwarg
         nside = 128
         nan_indx = tuple([0, 100])
-        scale = hp.nside2resol(nside) * 5
-        to_mask = utils._hp_grow_mask(nside, nan_indx, grow_size=scale, scale=None)
-
-        ra, dec = utils._hpid2_ra_dec(nside, np.arange(hp.nside2npix(nside)))
-        to_mask_set = set(to_mask)
-        all_close = []
-        for indx in nan_indx:
-            distances = utils._angular_separation(ra, dec, ra[indx], dec[indx])
-            close = np.where(distances <= scale)[0]
-            all_close.extend(close)
-            assert set(close).issubset(to_mask_set)
-        all_close = np.unique(all_close)
-        # Make sure we aren't including any extra pixels
-        assert np.size(all_close) == np.size(to_mask)
+        for scale in np.radians(
+            [
+                0.1,
+                1.0,
+                10.0,
+                30,
+            ]
+        ):
+            to_mask = utils._hp_grow_mask(nside, nan_indx, grow_size=scale, scale=None)
+
+            ra, dec = utils._hpid2_ra_dec(nside, np.arange(hp.nside2npix(nside)))
+            to_mask_set = set(to_mask)
+            all_close = []
+            for indx in nan_indx:
+                distances = utils._angular_separation(ra, dec, ra[indx], dec[indx])
+                close = np.where(distances <= scale)[0]
+                all_close.extend(close)
+                assert set(close).issubset(to_mask_set)
+            all_close = np.unique(all_close)
+            # Make sure we aren't including any extra pixels
+            assert np.size(all_close) == np.size(to_mask)
 
         # Check that 0 distance doesn't mask anything new
         to_mask = utils._hp_grow_mask(nside, nan_indx, grow_size=0)
         assert np.size(to_mask) == np.size(nan_indx)
+        # Check scale can be wiped out.
+        to_mask = utils._hp_grow_mask(nside, nan_indx, grow_size=0, scale=None)
+        assert np.size(to_mask) == np.size(nan_indx)
 
 
 if __name__ == "__main__":