# Licensed under a 3-clause BSD style license - see LICENSE.rst
"""
This module defines circular regions in both pixel and sky coordinates.
"""
import math
import astropy.units as u
import numpy as np
from astropy.coordinates import Angle
from regions._geometry import circular_overlap_grid
from regions._utils.spherical_helpers import (
get_circle_latitude_tangent_limits, get_circle_longitude_tangent_limits)
from regions._utils.wcs_helpers import (pixel_to_sky_mean_scale,
sky_to_pixel_mean_scale)
from regions.core.attributes import (PositiveScalar, PositiveScalarAngle,
RegionMetaDescr, RegionVisualDescr,
ScalarPixCoord, ScalarSkyCoord)
from regions.core.bounding_box import RegionBoundingBox
from regions.core.core import PixelRegion, SkyRegion, SphericalSkyRegion
from regions.core.mask import RegionMask
from regions.core.metadata import RegionMeta, RegionVisual
from regions.core.pixcoord import PixCoord
from regions.shapes.polygon import PolygonPixelRegion
__all__ = ['CirclePixelRegion', 'CircleSkyRegion', 'CircleSphericalSkyRegion']
[docs]
class CirclePixelRegion(PixelRegion):
"""
A circle defined using pixel coordinates.
Parameters
----------
center : `~regions.PixCoord`
The center position.
radius : float
The radius in pixels.
meta : `~regions.RegionMeta` or `dict`, optional
A dictionary that stores the meta attributes of the region.
visual : `~regions.RegionVisual` or `dict`, optional
A dictionary that stores the visual meta attributes of the
region.
Examples
--------
.. plot::
:include-source:
from regions import PixCoord, CirclePixelRegion
import matplotlib.pyplot as plt
fig, ax = plt.subplots(1, 1)
reg = CirclePixelRegion(PixCoord(x=8, y=7), radius=3.5)
patch = reg.plot(ax=ax, facecolor='none', edgecolor='red', lw=2,
label='Circle')
ax.legend(handles=(patch,), loc='upper center')
ax.set_xlim(0, 15)
ax.set_ylim(0, 15)
ax.set_aspect('equal')
"""
_params = ('center', 'radius')
_mpl_artist = 'Patch'
center = ScalarPixCoord('The center pixel position as a |PixCoord|.')
radius = PositiveScalar('The radius in pixels as a float.')
meta = RegionMetaDescr('The meta attributes as a |RegionMeta|')
visual = RegionVisualDescr('The visual attributes as a |RegionVisual|.')
def __init__(self, center, radius, meta=None, visual=None):
self.center = center
self.radius = radius
self.meta = meta or RegionMeta()
self.visual = visual or RegionVisual()
@property
def area(self):
return math.pi * self.radius ** 2
[docs]
def contains(self, pixcoord):
pixcoord = PixCoord._validate(pixcoord, name='pixcoord')
in_circle = self.center.separation(pixcoord) < self.radius
if self.meta.get('include', True):
return in_circle
else:
return np.logical_not(in_circle)
[docs]
def to_sky(self, wcs):
center, mean_scale = pixel_to_sky_mean_scale(self.center, wcs)
radius = Angle(self.radius * mean_scale, 'arcsec')
return CircleSkyRegion(center, radius, meta=self.meta.copy(),
visual=self.visual.copy())
[docs]
def to_spherical_sky(self, wcs=None, include_boundary_distortions=False,
n_points=None):
self._validate_planar_spherical_transform(wcs, include_boundary_distortions)
if include_boundary_distortions:
# Requires planar to spherical projection (using WCS) and discretization
# Will require implementing discretization in pixel space
# to get correct handling of distortions.
raise NotImplementedError
# ### Potential solution:
# # Leverage polygon class to_spherical_sky() functionality without
# # distortions, as the distortions were already computed in creating
# # that polygon approximation
# return self.discretize_boundary(n_points=npoints).to_spherical_sky(
# wcs=wcs, include_boundary_distortions=False
# )
return self.to_sky(wcs).to_spherical_sky()
@property
def bounding_box(self):
"""
Bounding box (`~regions.RegionBoundingBox`).
"""
xmin = self.center.x - self.radius
xmax = self.center.x + self.radius
ymin = self.center.y - self.radius
ymax = self.center.y + self.radius
return RegionBoundingBox.from_float(xmin, xmax, ymin, ymax)
[docs]
def to_mask(self, mode='center', subpixels=1):
self._validate_mode(mode, subpixels)
if mode == 'center':
mode = 'subpixels'
subpixels = 1
# Find bounding box and mask size
bbox = self.bounding_box
ny, nx = bbox.shape
# Find position of pixel edges and recenter so that circle is at
# origin
xmin = float(bbox.ixmin) - 0.5 - self.center.x
xmax = float(bbox.ixmax) - 0.5 - self.center.x
ymin = float(bbox.iymin) - 0.5 - self.center.y
ymax = float(bbox.iymax) - 0.5 - self.center.y
use_exact = 0 if mode == 'subpixels' else 1
fraction = circular_overlap_grid(xmin, xmax, ymin, ymax, nx, ny,
self.radius, use_exact, subpixels)
return RegionMask(fraction, bbox=bbox)
[docs]
def as_artist(self, origin=(0, 0), **kwargs):
"""
Return a matplotlib patch object for this region
(`matplotlib.patches.Circle`).
Parameters
----------
origin : array_like, optional
The ``(x, y)`` pixel position of the origin of the displayed
image.
**kwargs : dict
Any keyword arguments accepted by
`~matplotlib.patches.Circle`. These keywords will override
any visual meta attributes of this region.
Returns
-------
artist : `~matplotlib.patches.Circle`
A matplotlib circle patch.
"""
from matplotlib.patches import Circle
xy = self.center.x - origin[0], self.center.y - origin[1]
radius = self.radius
mpl_kwargs = self.visual.define_mpl_kwargs(self._mpl_artist)
mpl_kwargs.update(kwargs)
return Circle(xy=xy, radius=radius, **mpl_kwargs)
[docs]
def rotate(self, center, angle):
"""
Rotate the region.
Positive ``angle`` corresponds to counter-clockwise rotation.
Parameters
----------
center : `~regions.PixCoord`
The rotation center point.
angle : `~astropy.coordinates.Angle`
The rotation angle.
Returns
-------
region : `CirclePixelRegion`
The rotated region (which is an independent copy).
"""
center = self.center.rotate(center, angle)
return self.copy(center=center)
[docs]
def to_polygon(self, n_points=100):
"""
Return a `~regions.PolygonPixelRegion` that approximates this
circle.
Parameters
----------
n_points : int, optional
The number of polygon vertices. Default is 100.
Returns
-------
polygon : `~regions.PolygonPixelRegion`
A polygon region approximating the circle.
"""
theta = np.linspace(0, 2 * np.pi, n_points, endpoint=False)
x = self.radius * np.cos(theta) + self.center.x
y = self.radius * np.sin(theta) + self.center.y
vertices = PixCoord(x=x, y=y)
return PolygonPixelRegion(vertices=vertices,
meta=self.meta.copy(),
visual=self.visual.copy())
[docs]
class CircleSkyRegion(SkyRegion):
"""
A circle defined using sky coordinates.
Parameters
----------
center : `~astropy.coordinates.SkyCoord`
The center position.
radius : `~astropy.units.Quantity`
The radius in angular units.
meta : `~regions.RegionMeta` or `dict`, optional
A dictionary that stores the meta attributes of the region.
visual : `~regions.RegionVisual` or `dict`, optional
A dictionary that stores the visual meta attributes of the
region.
"""
_params = ('center', 'radius')
center = ScalarSkyCoord('The center position as a |SkyCoord|.')
radius = PositiveScalarAngle('The radius as a |Quantity| angle.')
meta = RegionMetaDescr('The meta attributes as a |RegionMeta|')
visual = RegionVisualDescr('The visual attributes as a |RegionVisual|.')
def __init__(self, center, radius, meta=None, visual=None):
self.center = center
self.radius = radius
self.meta = meta or RegionMeta()
self.visual = visual or RegionVisual()
[docs]
def to_polygon(self, wcs, n_points=100):
"""
Return a `~regions.PolygonSkyRegion` that approximates this
circle.
Parameters
----------
wcs : `~astropy.wcs.WCS`
The WCS to use for the sky-to-pixel-to-sky conversion.
n_points : int, optional
The number of polygon vertices. Default is 100.
Returns
-------
polygon : `~regions.PolygonSkyRegion`
A polygon region approximating the circle.
"""
return self.to_pixel(wcs).to_polygon(n_points=n_points).to_sky(wcs)
[docs]
def to_pixel(self, wcs):
center, mean_scale = sky_to_pixel_mean_scale(self.center, wcs)
radius = self.radius.to(u.arcsec).value * mean_scale
return CirclePixelRegion(center, radius, meta=self.meta.copy(),
visual=self.visual.copy())
[docs]
def to_spherical_sky(self, wcs=None, include_boundary_distortions=False,
n_points=None):
self._validate_planar_spherical_transform(wcs, include_boundary_distortions)
if include_boundary_distortions:
# Requires planar to spherical projection (using WCS) and discretization
# Will require implementing discretization in pixel space
# to get correct handling of distortions.
raise NotImplementedError
# ### Potential solution:
# # Leverage polygon class to_spherical_sky() functionality without
# # distortions, as the distortions were already computed in creating
# # that polygon approximation
# return self.to_pixel(wcs).discretize_boundary(n_points=n_points).to_spherical_sky(
# wcs=wcs, include_boundary_distortions=False
# )
return CircleSphericalSkyRegion(
self.center, self.radius, meta=self.meta, visual=self.visual,
)
[docs]
class CircleSphericalSkyRegion(SphericalSkyRegion):
"""
Class for a circular sky region, where the circle is interpreted
within a spherical geometry reference frame.
Parameters
----------
center : `~astropy.coordinates.SkyCoord`
The center position.
radius : `~astropy.units.Quantity`
The radius in angular units.
meta : `~regions.RegionMeta` or `dict`, optional
A dictionary that stores the meta attributes of the region.
visual : `~regions.RegionVisual` or `dict`, optional
A dictionary that stores the visual meta attributes of the
region.
"""
_params = ('center', 'radius')
center = ScalarSkyCoord('The center position as a |SkyCoord|.')
radius = PositiveScalarAngle('The radius as a |Quantity| angle.')
meta = RegionMetaDescr('The meta attributes as a |RegionMeta|')
visual = RegionVisualDescr('The visual attributes as a |RegionVisual|.')
def __init__(self, center, radius, meta=None, visual=None):
self.center = center
self.radius = radius
self.meta = meta or RegionMeta()
self.visual = visual or RegionVisual()
[docs]
def contains(self, coord):
in_circle = self.center.separation(coord) < self.radius
if self.meta.get('include', True):
return in_circle
else:
return np.logical_not(in_circle)
@property
def bounding_circle(self):
return self.copy()
@property
def bounding_lonlat(self):
lons_arr = get_circle_longitude_tangent_limits(self.center, self.radius)
lats_arr = get_circle_latitude_tangent_limits(self.center, self.radius)
lons_arr, lats_arr = self._validate_lonlat_bounds(lons_arr, lats_arr)
return lons_arr, lats_arr
[docs]
def discretize_boundary(self, n_points=100):
# Avoid circular imports:
from .polygon import PolygonSphericalSkyRegion
theta = np.linspace(0, 1, num=n_points, endpoint=False) * 360 * u.deg
# Need to invert order because of CW convention:
bound_verts = self.center.directional_offset_by(theta[::-1], self.radius)
return PolygonSphericalSkyRegion(bound_verts, meta=self.meta.copy(),
visual=self.visual.copy())
[docs]
def to_polygon(self, n_points=100):
"""
Return a `~regions.PolygonSphericalSkyRegion` that approximates this
circle.
Parameters
----------
n_points : int, optional
The number of polygon vertices. Default is 100.
Returns
-------
polygon : `~regions.PolygonSphericalSkyRegion`
A polygon region approximating the circle.
"""
return self.discretize_boundary(n_points=n_points)
[docs]
def to_sky(
self, wcs=None, include_boundary_distortions=False, n_points=None,
):
self._validate_planar_spherical_transform(wcs, include_boundary_distortions)
if include_boundary_distortions:
# Requires spherical to planar projection (from WCS) and discretization
# Use to_pixel(), then apply "small angle approx" to get planar sky.
return self.to_pixel(
include_boundary_distortions=include_boundary_distortions,
wcs=wcs,
n_points=n_points,
).to_sky(wcs)
return CircleSkyRegion(
self.center, self.radius, meta=self.meta, visual=self.visual,
)
[docs]
def to_pixel(
self,
wcs=None,
include_boundary_distortions=False,
n_points=None,
):
self._validate_planar_spherical_transform(wcs, include_boundary_distortions)
if include_boundary_distortions:
from .polygon import PolygonPixelRegion
disc_kwargs = {} if n_points is None else {'n_points': n_points}
# Requires spherical to planar projection (from WCS) and discretization
verts = wcs.world_to_pixel(
self.discretize_boundary(**disc_kwargs).vertices,
)
return PolygonPixelRegion(
PixCoord(*verts), meta=self.meta.copy(), visual=self.visual.copy(),
)
return self.to_sky().to_pixel(wcs)