Introduce RegionGeom.from_regions()
@@ -6,11 +6,12 @@ from astropy.coordinates import SkyCoord
from astropy.io import fits
from astropy.table import Table
from astropy.wcs.utils import proj_plane_pixel_area, wcs_to_celestial_frame
-from regions import FITSRegionParser, fits_region_objects_to_table
+from regions import FITSRegionParser, fits_region_objects_to_table, SkyRegion, CompoundSkyRegion
from gammapy.utils.regions import (
compound_region_to_list,
list_to_compound_region,
make_region,
+ compound_region_center
)
from gammapy.maps.wcs import _check_width
from .core import MapCoord, Map
@@ -53,11 +54,16 @@ class RegionGeom(Geom):
self._axes = MapAxes.from_default(axes)
if wcs is None and region is not None:
+ if isinstance(region, CompoundSkyRegion):
+ center = compound_region_center(region)
+ else:
+ center = region.center
+
wcs = WcsGeom.create(
- skydir=region.center,
+ skydir=center,
binsz=self.binsz,
proj=self.projection,
- frame=self.frame,
+ frame=center.frame.name,
).wcs
self._wcs = wcs
@@ -519,6 +525,28 @@ class RegionGeom(Geom):
return hdulist
+ @classmethod
+ def from_regions(cls, regions, **kwargs):
+ """Create region geom from list of regions
+
+ The regions are combined to a compound region.
+
+ Parameters
+ ----------
+ regions : list of `~regions.SkyRegion`
+ Regions
+
+ Returns
+ -------
+ geom : `RegionGeom`
+ Region map geometry
+ """
+ if isinstance(regions, SkyRegion):
+ regions = [regions]
+
+ region = list_to_compound_region(regions)
+ return cls(region, **kwargs)
+
@classmethod
def from_hdulist(cls, hdulist, format="ogip"):
"""Read region table and convert it to region list.
