Sky Map Plotting (ligo.skymap.plot.allsky)

Axes subclasses for all-sky maps. This adds several astropy.visualization.wcsaxes.WCSAxes subclasses to the Matplotlib projection registry. The projections are:

  • astro degrees aitoff

  • astro degrees mollweide

  • astro hours aitoff

  • astro hours mollweide

  • geo degrees aitoff

  • geo hours aitoff

  • geo degrees mollweide

  • geo hours mollweide

  • astro globe with options center and rotate

  • astro zoom with options center, radius, and rotate

Example

The following example demonstrates most of the features of this module.

from astropy.coordinates import SkyCoord
from astropy.io import fits
from astropy import units as u
import ligo.skymap.plot
from matplotlib import pyplot as plt

url = 'https://dcc.ligo.org/public/0146/G1701985/001/bayestar_no_virgo.fits.gz'
center = SkyCoord.from_name('NGC 4993')

fig = plt.figure(figsize=(4, 4), dpi=100)

ax = plt.axes(
    [0.05, 0.05, 0.9, 0.9],
    projection='astro globe',
    center=center)

ax_inset = plt.axes(
    [0.59, 0.3, 0.4, 0.4],
    projection='astro zoom',
    center=center,
    radius=10*u.deg)

for key in ['ra', 'dec']:
    ax_inset.coords[key].set_ticklabel_visible(False)
    ax_inset.coords[key].set_ticks_visible(False)
ax.grid()
ax.mark_inset_axes(ax_inset)
ax.connect_inset_axes(ax_inset, 'upper left')
ax.connect_inset_axes(ax_inset, 'lower left')
ax_inset.scalebar((0.1, 0.1), 5 * u.deg).label()
ax_inset.compass(0.9, 0.1, 0.2)

ax.imshow_hpx(url, cmap='cylon')
ax_inset.imshow_hpx(url, cmap='cylon')
ax_inset.plot(
    center.ra.deg, center.dec.deg,
    transform=ax_inset.get_transform('world'),
    marker=ligo.skymap.plot.reticle(),
    markersize=30,
    markeredgewidth=3)

(Source code, png, hires.png, pdf)

../../../_images/allsky-1.png
class ligo.skymap.plot.allsky.AutoScaledWCSAxes(*args, header, **kwargs)[source]

Bases: astropy.visualization.wcsaxes.WCSAxes

Axes base class. The pixel scale is adjusted to the DPI of the image, and there are a variety of convenience methods.

Build an axes in a figure.

Parameters
figFigure

The axes is build in the Figure fig.

rect[left, bottom, width, height]

The axes is build in the rectangle rect. rect is in Figure coordinates.

sharex, shareyAxes, optional

The x or y axis is shared with the x or y axis in the input Axes.

frameonbool, optional

True means that the axes frame is visible.

**kwargs
Other optional keyword arguments:

adjustable: {‘box’, ‘datalim’}

agg_filter: a filter function, which takes a (m, n, 3) float array and a dpi value, and returns a (m, n, 3) array alpha: float anchor: 2-tuple of floats or {‘C’, ‘SW’, ‘S’, ‘SE’, …} animated: bool aspect: {‘auto’, ‘equal’} or num autoscale_on: bool autoscalex_on: bool autoscaley_on: bool axes_locator: Callable[[Axes, Renderer], Bbox] axisbelow: bool or ‘line’ clip_box: Bbox clip_on: bool clip_path: [(Path, Transform) | Patch | None] contains: callable facecolor: color fc: color figure: Figure frame_on: bool gid: str in_layout: bool label: object navigate: bool navigate_mode: unknown path_effects: AbstractPathEffect picker: None or bool or float or callable position: [left, bottom, width, height] or Bbox rasterization_zorder: float or None rasterized: bool or None sketch_params: (scale: float, length: float, randomness: float) snap: bool or None title: str transform: Transform url: str visible: bool xbound: unknown xlabel: str xlim: (left: float, right: float) xmargin: float greater than -0.5 xscale: {“linear”, “log”, “symlog”, “logit”, …} xticklabels: List[str] xticks: list ybound: unknown ylabel: str ylim: (bottom: float, top: float) ymargin: float greater than -0.5 yscale: {“linear”, “log”, “symlog”, “logit”, …} yticklabels: List[str] yticks: list zorder: float

Returns
axesAxes

The new Axes object.

compass(self, x, y, size)[source]

Add a compass to indicate the north and east directions.

Parameters
x, yfloat

Position of compass vertex in axes coordinates.

sizefloat

Size of compass in axes coordinates.

connect_inset_axes(self, ax, loc, *args, **kwargs)[source]

Convenience function to connect a corner of another WCSAxes to the matching point inside this one.

Parameters
axastropy.visualization.wcsaxes.WCSAxes

The other axes.

locint, str

Which corner to connect. For valid values, see matplotlib.offsetbox.AnchoredOffsetbox.

Returns
patchmatplotlib.patches.ConnectionPatch
Other Parameters
args :

Extra arguments for matplotlib.patches.ConnectionPatch

kwargs :

Extra keyword arguments for matplotlib.patches.ConnectionPatch

contour_hpx(self, data, hdu_in=None, order='bilinear', nested=False, field=0, smooth=None, **kwargs)[source]

Add contour levels for a HEALPix data set.

Parameters
datanumpy.ndarray or str or TableHDU or BinTableHDU or tuple

The HEALPix data set. If this is a numpy.ndarray, then it is interpreted as the HEALPi array in the same coordinate system as the axes. Otherwise, the input data can be any type that is understood by reproject.reproject_from_healpix.

smoothastropy.units.Quantity, optional

An optional smoothing length in angle-compatible units.

Returns
countoursmatplotlib.contour.QuadContourSet
Other Parameters
hdu_in, order, nested, field, smooth :

Extra arguments for reproject.reproject_from_healpix

kwargs :

Extra keyword arguments for matplotlib.axes.Axes.contour

contourf_hpx(self, data, hdu_in=None, order='bilinear', nested=False, field=0, smooth=None, **kwargs)[source]

Add filled contour levels for a HEALPix data set.

Parameters
datanumpy.ndarray or str or TableHDU or BinTableHDU or tuple

The HEALPix data set. If this is a numpy.ndarray, then it is interpreted as the HEALPi array in the same coordinate system as the axes. Otherwise, the input data can be any type that is understood by reproject.reproject_from_healpix.

smoothastropy.units.Quantity, optional

An optional smoothing length in angle-compatible units.

Returns
contoursmatplotlib.contour.QuadContourSet
Other Parameters
hdu_in, order, nested, field, smooth :

Extra arguments for reproject.reproject_from_healpix

kwargs :

Extra keyword arguments for matplotlib.axes.Axes.contour

imshow_hpx(self, data, hdu_in=None, order='bilinear', nested=False, field=0, smooth=None, **kwargs)[source]

Add an image for a HEALPix data set.

Parameters
datanumpy.ndarray or str or TableHDU or BinTableHDU or tuple

The HEALPix data set. If this is a numpy.ndarray, then it is interpreted as the HEALPi array in the same coordinate system as the axes. Otherwise, the input data can be any type that is understood by reproject.reproject_from_healpix.

smoothastropy.units.Quantity, optional

An optional smoothing length in angle-compatible units.

Returns
imagematplotlib.image.AxesImage
Other Parameters
hdu_in, order, nested, field, smooth :

Extra arguments for reproject.reproject_from_healpix

kwargs :

Extra keyword arguments for matplotlib.axes.Axes.contour

mark_inset_axes(self, ax, *args, **kwargs)[source]

Outline the footprint of another WCSAxes inside this one.

Parameters
axastropy.visualization.wcsaxes.WCSAxes

The other axes.

Returns
patchmatplotlib.patches.PathPatch
Other Parameters
args :

Extra arguments for matplotlib.patches.PathPatch

kwargs :

Extra keyword arguments for matplotlib.patches.PathPatch

scalebar(self, *args, **kwargs)[source]

Add scale bar.

Parameters
xytuple

The axes coordinates of the scale bar.

lengthastropy.units.Quantity

The length of the scale bar in angle-compatible units.

Returns
patchmatplotlib.patches.FancyArrowPatch
Other Parameters
args :

Extra arguments for matplotlib.patches.FancyArrowPatch

kwargs :

Extra keyword arguments for matplotlib.patches.FancyArrowPatch