Merge 9005878 into 7116cd8

hips/draw/simple.py

@@ -1,16 +1,15 @@
 # Licensed under a 3-clause BSD style license - see LICENSE.rst
 """HiPS tile drawing -- simple method."""
-from typing import Tuple, List
 import numpy as np
+from typing import Tuple, List
 from astropy.wcs.utils import proj_plane_pixel_scales
 from skimage.transform import ProjectiveTransform, warp
 from ..tiles import HipsSurveyProperties, HipsTile, HipsTileMeta
 from ..utils import WCSGeometry, compute_healpix_pixel_indices, get_hips_order_for_resolution
 
 __all__ = [
     'make_sky_image',
-    'SimpleTilePainter',
-
+    'SimpleTilePainter'
 ]
 
 __doctest_skip__ = [
@@ -178,11 +177,49 @@ def draw_tiles(self) -> np.ndarray:
 
         return image
 
+    def draw_hips_tile_grid(self) -> None:
+        """Draw lines on the output image (mainly used for debugging)."""
+        import matplotlib.pyplot as plt
+        for tile in self.tiles:
+            corners = tile.meta.skycoord_corners.transform_to(self.geometry.celestial_frame)
+            ax = plt.subplot(projection=self.geometry.wcs)
+            ax.plot(corners.data.lon.deg, corners.data.lat.deg,
+                    'red', lw=1, transform=ax.get_transform('icrs'))
+        ax.imshow(self.image, origin='lower')
+
     def run(self) -> None:
         """Run all steps of the naive algorithm."""
         self.float_image = self.draw_tiles()
 
 
+def draw_debug_image(geometry: WCSGeometry, tile: HipsTile, image: np.ndarray) -> None:
+    """Draw markers on the output image (mainly used for debugging).
+
+    The following denotes their correspondence:
+    * red <=> North
+    * green <=> West
+    * blue <=> South
+    * yellow <=> East
+
+    Parameters
+    ----------
+    geometry : `~hips.utils.WCSGeometry`
+        Geometry of the output image
+    tile : HipsTile
+        HiPS tile
+    image : np.ndarray
+        Image containing HiPS tiles
+    """
+    import matplotlib.pyplot as plt
+    corners = tile.meta.skycoord_corners.transform_to(geometry.celestial_frame)
+    colors = ['red', 'green', 'blue', 'yellow']
+    ax = plt.subplot(projection=geometry.wcs)
+    for index, corner in enumerate(corners):
+        ax.scatter(corner.data.lon.deg, corner.data.lat.deg,
+                   s=80, transform=ax.get_transform('icrs'), color=colors[index])
+    ax.imshow(image, origin='lower')
+
+
 def make_sky_image(geometry: WCSGeometry, hips_survey: HipsSurveyProperties, tile_format: str) -> np.ndarray:
     """Make sky image: fetch tiles and draw.
 

hips/draw/tests/test_simple.py

@@ -5,7 +5,7 @@
 from astropy.coordinates import SkyCoord
 from astropy.tests.helper import remote_data
 from ...tiles import HipsSurveyProperties
-from ..simple import make_sky_image, SimpleTilePainter
+from ..simple import make_sky_image, SimpleTilePainter, draw_debug_image
 from ...utils.wcs import WCSGeometry
 from ...utils.testing import make_test_wcs_geometry, requires_hips_extra
 
@@ -95,3 +95,9 @@ def test_compute_matching_hips_order(self, pars):
     def test_run(self):
         self.simple_tile_painter.run()
         assert_allclose(self.simple_tile_painter.image[200, 994], 2120)
+
+    def test_draw_hips_tile_grid(self):
+        self.simple_tile_painter.draw_hips_tile_grid()
+
+    def test_draw_debug_image(self):
+        draw_debug_image(self.geometry, self.simple_tile_painter.tiles[3], self.simple_tile_painter.image)