Remove SkyMask (merge with SkyImage)

docs/background/make_reflected_regions.py

@@ -1,26 +1,29 @@
 """Example how to compute and plot reflected regions."""
 from astropy.coordinates import SkyCoord, Angle
 from regions import CircleSkyRegion
-from gammapy.image import SkyMask
+from gammapy.image import SkyImage
 from gammapy.background import find_reflected_regions
 
-mask = SkyMask.empty(name='Exclusion Mask', nxpix=801, nypix=701, binsz=0.01,
-                     coordsys='CEL', xref=83.2, yref=22.7)
-mask.fill_random_circles(n=8, min_rad=30, max_rad=80)
+mask = SkyImage.empty(
+    name='Exclusion Mask', nxpix=801, nypix=701, binsz=0.01,
+    coordsys='CEL', xref=83.633, yref=23.014, fill=1,
+)
 
-pos = SkyCoord(80.2, 23.5, unit='deg')
-radius = Angle(0.4, 'deg')
-test_region = CircleSkyRegion(pos, radius)
-center = SkyCoord(82.8, 22.5, unit='deg')
-regions = find_reflected_regions(test_region, center, mask)
+pos = SkyCoord(83.633, 22.014, unit='deg')
+radius = Angle(0.3, 'deg')
+on_region = CircleSkyRegion(pos, radius)
+center = SkyCoord(83.633, 23.014, unit='deg')
+regions = find_reflected_regions(on_region, center, mask)
 
 import matplotlib.pyplot as plt
 
-fig = plt.figure(figsize=(8, 5), dpi=80)
+fig = plt.figure(figsize=(6, 5))
 ax = fig.add_axes([0.1, 0.1, 0.8, 0.8], projection=mask.wcs)
 mask.plot(ax, fig)
 for reg in regions:
-    reg.to_pixel(mask.wcs).plot(ax, facecolor='red')
-test_region.to_pixel(mask.wcs).plot(ax, facecolor='blue')
+    reg.to_pixel(mask.wcs).plot(ax, color='r')
+on_region.to_pixel(mask.wcs).plot(ax)
+ax.scatter(center.ra.deg, center.dec.deg, s=50, c='w', marker='+',
+           transform=ax.get_transform('world'))
 
 plt.show()

examples/test_image_analysis.py

@@ -4,7 +4,7 @@
 from astropy.coordinates import SkyCoord, Angle
 from gammapy.utils.energy import Energy
 from gammapy.data import DataStore
-from gammapy.image import SkyImage, SkyMask
+from gammapy.image import SkyImage
 from gammapy.background import OffDataBackgroundMaker
 from gammapy.scripts import StackedObsImageMaker
 import matplotlib.pyplot as plt
@@ -63,7 +63,7 @@ def make_image_from_2d_bg():
                            yref=center.b.deg, proj='TAN', coordsys='GAL')
 
     refheader = image.to_image_hdu().header
-    exclusion_mask = SkyMask.read('$GAMMAPY_EXTRA/datasets/exclusion_masks/tevcat_exclusion.fits')
+    exclusion_mask = SkyImage.read('$GAMMAPY_EXTRA/datasets/exclusion_masks/tevcat_exclusion.fits')
     exclusion_mask = exclusion_mask.reproject(reference=refheader)
     mosaic_images = StackedObsImageMaker(image, energy_band=energy_band, offset_band=offset_band, data_store=data_store,
                                          obs_table=data_store.obs_table, exclusion_mask=exclusion_mask)

gammapy/background/reflected.py

@@ -3,21 +3,19 @@
 import numpy as np
 from astropy.coordinates import Angle
 from regions import PixCoord, CirclePixelRegion
-from ..image import SkyMask
+from ..image import SkyImage
 from .background_estimate import BackgroundEstimate
-import logging
 
 __all__ = [
     'find_reflected_regions',
     'ReflectedRegionsBackgroundEstimator',
 ]
 
-log = logging.getLogger(__name__)
-
 
 def find_reflected_regions(region, center, exclusion_mask=None,
-                           angle_increment=None, min_distance=None,
-                           min_distance_input=None):
+                           angle_increment='0.1 rad',
+                           min_distance='0 rad',
+                           min_distance_input='0.1 rad'):
     """Find reflected regions.
 
     Converts to pixel coordinates internally.
@@ -28,36 +26,29 @@ def find_reflected_regions(region, center, exclusion_mask=None,
         Region
     center : `~astropy.coordinates.SkyCoord`
         Rotation point
-    exclusion_mask : `~gammapy.image.SkyMask`, optional
+    exclusion_mask : `~gammapy.image.SkyImage`, optional
         Exclusion mask
     angle_increment : `~astropy.coordinates.Angle`
-        Rotation angle for each step, default: 0.1 rad
+        Rotation angle for each step
     min_distance : `~astropy.coordinates.Angle`
-        Minimal distance between to reflected regions, default: 0 rad
+        Minimal distance between to reflected regions
     min_distance_input : `~astropy.coordinates.Angle`
-        Minimal distance from input region, default: 0.1 rad
+        Minimal distance from input region
 
     Returns
     -------
     regions : list of `~regions.SkyRegion`
         Reflected regions list
     """
-    angle_increment = Angle('0.1 rad') if angle_increment is None else Angle(angle_increment)
-    min_distance = Angle('0 rad') if min_distance is None else Angle(min_distance)
-    min_distance_input = Angle('0.1 rad') if min_distance_input is None else Angle(min_distance_input)
+    angle_increment = Angle(angle_increment)
+    min_distance = Angle(min_distance)
+    min_distance_input = Angle(min_distance_input)
 
     # Create empty exclusion mask if None is provided
     if exclusion_mask is None:
-        min_size = region.center.separation(center)
-        binsz = 0.02
-        npix = int((3 * min_size / binsz).value)
-        exclusion_mask = SkyMask.empty(name='empty exclusion mask',
-                                       xref=center.galactic.l.value,
-                                       yref=center.galactic.b.value,
-                                       binsz=binsz,
-                                       nxpix=npix,
-                                       nypix=npix,
-                                       fill=1)
+        exclusion_mask = make_default_exclusion_mask(center, exclusion_mask, region)
+
+    distance_image = exclusion_mask.distance_image
 
     reflected_regions_pix = list()
     wcs = exclusion_mask.wcs
@@ -77,26 +68,40 @@ def find_reflected_regions(region, center, exclusion_mask=None,
     # Add required minimal distance between two off regions
     min_ang += min_distance
 
-    # Maximum allowd angle before the an overlap with the ON regions happens
-    max_angle = angle + Angle('360deg') - min_ang - min_distance_input
+    # Maximum allowed angle before the an overlap with the ON regions happens
+    max_angle = angle + Angle('360 deg') - min_ang - min_distance_input
 
     # Starting angle
     curr_angle = angle + min_ang + min_distance_input
+
     while curr_angle < max_angle:
         test_pos = _compute_xy(pix_center, offset, curr_angle)
         test_reg = CirclePixelRegion(test_pos, pix_region.radius)
-        if not _is_inside_exclusion(test_reg, exclusion_mask):
+        if distance_image.lookup_pix(test_reg.center) > pix_region.radius:
             reflected_regions_pix.append(test_reg)
             curr_angle = curr_angle + min_ang
         else:
             curr_angle = curr_angle + angle_increment
 
     reflected_regions = [_.to_sky(wcs) for _ in reflected_regions_pix]
-    log.debug('Found {} reflected regions:\n {}'.format(len(reflected_regions),
-                                                        reflected_regions))
     return reflected_regions
 
 
+def make_default_exclusion_mask(center, region):
+    min_size = region.center.separation(center)
+    binsz = 0.02
+    npix = int((3 * min_size / binsz).value)
+    return SkyImage.empty(
+        name='empty exclusion mask',
+        xref=center.galactic.l.value,
+        yref=center.galactic.b.value,
+        binsz=binsz,
+        nxpix=npix,
+        nypix=npix,
+        fill=1,
+    )
+
+
 def _compute_xy(pix_center, offset, angle):
     """Compute x, y position for a given position angle and offset.
 
@@ -109,15 +114,6 @@ def _compute_xy(pix_center, offset, angle):
     return PixCoord(x=x, y=y)
 
 
-# TODO :Copied from gammapy.region.PixCircleList (deleted), find better place
-def _is_inside_exclusion(pixreg, exclusion):
-    x, y = pixreg.center.x, pixreg.center.y
-    image = exclusion.distance_image
-    excl_dist = image.data
-    val = excl_dist[np.round(y).astype(int), np.round(x).astype(int)]
-    return val < pixreg.radius
-
-
 class ReflectedRegionsBackgroundEstimator(object):
     """Reflected Regions background estimator.
 
@@ -127,7 +123,7 @@ class ReflectedRegionsBackgroundEstimator(object):
         Target region
     obs_list : `~gammapy.data.ObservationList`
         List of observations to process
-    exclusion : `~gammapy.image.SkyMask`
+    exclusion : `~gammapy.image.SkyImage`
         Exclusion mask
     config : dict
         Config dict to be passed to :func:`gammapy.background.find_reflected_regions`

gammapy/background/tests/test_reflected.py

@@ -1,19 +1,19 @@
 # Licensed under a 3-clause BSD style license - see LICENSE.rst
 from __future__ import absolute_import, division, print_function, unicode_literals
-from astropy.tests.helper import pytest
-from astropy.tests.helper import assert_quantity_allclose
+from astropy.tests.helper import pytest, assert_quantity_allclose
 from astropy.coordinates import SkyCoord, Angle
 from regions import CircleSkyRegion
-from ..reflected import find_reflected_regions, ReflectedRegionsBackgroundEstimator
-from ...image import SkyMask
 from ...utils.testing import requires_data, requires_dependency
-from ...data import Target, DataStore
+from ...image import SkyImage
+from ...data import DataStore
+from ..reflected import find_reflected_regions, ReflectedRegionsBackgroundEstimator
+
 
 @pytest.fixture
 def mask():
     """Example mask for testing."""
     filename = '$GAMMAPY_EXTRA/datasets/exclusion_masks/tevcat_exclusion.fits'
-    return SkyMask.read(filename, hdu='EXCLUSION')
+    return SkyImage.read(filename, hdu='EXCLUSION')
 
 @pytest.fixture
 def on_region():
@@ -23,6 +23,7 @@ def on_region():
     region = CircleSkyRegion(pos, radius)
     return region
 
+
 @pytest.fixture
 def obs_list():
     """Example observation list for testing."""
@@ -31,10 +32,11 @@ def obs_list():
     obs_ids = [23523, 23526]
     return datastore.obs_list(obs_ids)
 
+
 @requires_dependency('scipy')
 @requires_data('gammapy-extra')
 def test_find_reflected_regions(mask, on_region):
-    pointing = SkyCoord(83.2, 22.7, unit='deg', frame='icrs')
+    pointing = SkyCoord(83.2, 22.7, unit='deg')
     regions = find_reflected_regions(region=on_region, center=pointing,
                                      exclusion_mask=mask,
                                      min_distance_input=Angle('0 deg'))
@@ -46,9 +48,9 @@ def test_find_reflected_regions(mask, on_region):
 @requires_dependency('scipy')
 class TestReflectedRegionBackgroundEstimator:
     def setup(self):
-        temp = ReflectedRegionsBackgroundEstimator(on_region = on_region(),
-                                                   exclusion = mask(),
-                                                   obs_list = obs_list())
+        temp = ReflectedRegionsBackgroundEstimator(on_region=on_region(),
+                                                   exclusion=mask(),
+                                                   obs_list=obs_list())
         self.bg_maker = temp
 
     def test_basic(self):
@@ -61,7 +63,7 @@ def test_process(self, on_region, obs_list, mask):
         assert len(bg_estimate.off_region) == 22
 
     def test_run(self):
-        self.bg_maker.config.update(min_distance = '0.2 deg')
+        self.bg_maker.config.update(min_distance='0.2 deg')
         self.bg_maker.run()
         assert len(self.bg_maker.result[1].off_region) == 22
 
@@ -71,4 +73,3 @@ def test_plot(self):
         self.bg_maker.plot()
         self.bg_maker.plot(idx=1)
         self.bg_maker.plot(idx=[0, 1])
-

gammapy/cube/cube_pipe.py

@@ -40,7 +40,7 @@ class SingleObsCubeMaker(object):
     offset_band : `astropy.coordinates.Angle`
         Offset band selection
     exclusion_mask : `~gammapy.cube.SkyCube`
-        Cube of `~gammapy.image.SkyMask`
+        Exclusion mask
     save_bkg_scale: bool
         True if you want to save the normalisation of the bkg computed outside the exclusion region in a Table
     """
@@ -184,7 +184,7 @@ class StackedObsCubeMaker(object):
     obs_table : `~astropy.table.Table`
         Required columns: OBS_ID
     exclusion_mask : `~gammapy.cube.SkyCube`
-        Cube of `~gammapy.image.SkyMask`
+        Exclusion mask
     save_bkg_scale: bool
         True if you want to save the normalisation of the bkg for each run in a `Table` table_bkg_norm with two columns:
          "OBS_ID" and "bkg_scale"

gammapy/cube/tests/test_cube_pipe.py

@@ -1,25 +1,24 @@
 """Example how to make a Cube analysis from a 2D background model.
 """
-
 from __future__ import absolute_import, division, print_function, unicode_literals
 import numpy as np
 from numpy.testing import assert_allclose
 from astropy.tests.helper import assert_quantity_allclose
 from astropy.coordinates import SkyCoord, Angle
 from astropy.units import Quantity
-from gammapy.extern.pathlib import Path
-from ...data import DataStore
+from ...extern.pathlib import Path
 from ...utils.testing import requires_dependency, requires_data, pytest
-from ...image import SkyMask
+from ...utils.energy import Energy
+from ...data import DataStore
+from ...image import SkyImage
 from .. import StackedObsCubeMaker
 from ...background import OffDataBackgroundMaker
 from .. import SkyCube
-from ...utils.energy import Energy
 
 
 def make_empty_cube(image_size, energy, center, data_unit=None):
     """
-    Make an empty `SkyCube` from a given `SkyMask` and an energy binning.
+    Make an empty `SkyCube` from a given `SkyImage` and an energy binning.
 
     Parameters
     ----------
@@ -90,7 +89,7 @@ def test_cube_pipe(tmpdir):
     data_store = DataStore.from_dir(tmpdir)
 
     refheader = ref_cube_images.sky_image_ref.to_image_hdu().header
-    exclusion_mask = SkyMask.read('$GAMMAPY_EXTRA/datasets/exclusion_masks/tevcat_exclusion.fits')
+    exclusion_mask = SkyImage.read('$GAMMAPY_EXTRA/datasets/exclusion_masks/tevcat_exclusion.fits')
     exclusion_mask = exclusion_mask.reproject(reference=refheader)
     ref_cube_skymask.data = np.tile(exclusion_mask.data, (5, 1, 1))
     # Pb with the load psftable for one of the run that is not implemented yet...

gammapy/data/tests/test_obs_stats.py

@@ -8,7 +8,7 @@
 from ...data import DataStore, ObservationList, ObservationStats, Target
 from ...utils.testing import requires_data, requires_dependency
 from ...background import ReflectedRegionsBackgroundEstimator
-from ...image import SkyMask
+from ...image import SkyImage
 
 
 def get_obs_list():
@@ -36,7 +36,7 @@ def target():
 
 @pytest.fixture(scope='session')
 def mask():
-    return SkyMask.read('$GAMMAPY_EXTRA/datasets/exclusion_masks/tevcat_exclusion.fits')
+    return SkyImage.read('$GAMMAPY_EXTRA/datasets/exclusion_masks/tevcat_exclusion.fits')
 
 
 @pytest.fixture(scope='session')

gammapy/data/tests/test_obs_summary.py

@@ -9,7 +9,7 @@
 from ...data import Target
 from ...utils.testing import requires_data, requires_dependency
 from ...background import ReflectedRegionsBackgroundEstimator
-from ...image import SkyMask
+from ...image import SkyImage
 
 
 def table_summary():
@@ -52,7 +52,7 @@ def obs_summary():
     target = Target(position=pos, on_region=on_region,
                     name='Crab Nebula', tag='crab')
 
-    mask = SkyMask.read('$GAMMAPY_EXTRA/datasets/exclusion_masks/tevcat_exclusion.fits')
+    mask = SkyImage.read('$GAMMAPY_EXTRA/datasets/exclusion_masks/tevcat_exclusion.fits')
 
     obs_list = ObservationList([datastore.obs(_) for _ in run_list])
     obs_stats = ObservationStatsList()