mocpy.MOC

class mocpy.MOC(interval_set=None)[source]

Multi-order spatial coverage class

A MOC describes the coverage of an arbitrary region on the unit sphere. MOCs are usually used for describing the global coverage of catalog/image surveys such as GALEX or SDSS. A MOC corresponds to a list of HEALPix cells at different depths. This class gives you the possibility to:

  1. Define MOC objects:
  • From a FITS file that stores HEALPix cells (see from_fits).
  • Directly from a list of HEALPix cells expressed either as a numpy structural array (see from_healpix_cells) or a simple python dictionnary (see from_json).
  • From a list of sky coordinates (see from_skycoords, from_lonlat).
  • From a convex/concave polygon (see from_polygon).
  • From a cone (will be implemented in a next version).
  1. Perform fast logical operations between MOC objects:
  1. Plot the MOC objects:
  • Draw the MOC with its HEALPix cells (see fill)
  • Draw the perimeter of a MOC (see border)
  1. Get the sky coordinates defining the border(s) of MOC objects (see get_boundaries).
  2. Serialize MOC objects to astropy.io.fits.HDUList or JSON dictionary and save it to a file.
Attributes:
LARK_PARSER_STR
max_order

Depth of the smallest HEALPix cells found in the MOC instance.

sky_fraction

Sky fraction covered by the MOC

Methods

add_neighbours(self) Extends the MOC instance so that it includes the HEALPix cells touching its border.
border(self, ax, wcs, \*\*kw_mpl_pathpatch) Draws the MOC border(s) on a matplotlib axis.
complement(self) Returns the complement of the MOC instance.
contains(self, ra, dec[, keep_inside]) Returns a boolean mask array of the positions lying inside (or outside) the MOC instance.
degrade_to_order(self, new_order) Degrades the MOC instance to a new, less precise, MOC.
difference(self, another_moc, \*args) Difference between the MOC instance and other MOCs.
empty(self) Checks whether the MOC is empty.
fill(self, ax, wcs, \*\*kw_mpl_pathpatch) Draws the MOC on a matplotlib axis.
from_fits(filename) Loads a MOC from a FITS file.
from_fits_images(path_l, max_norder) Loads a MOC from a set of FITS file images.
from_healpix_cells(ipix, depth[, fully_covered]) Creates a MOC from a set of HEALPix cells at a given depth.
from_image(header, max_norder[, mask]) Creates a MOC from an image stored as a FITS file.
from_ivorn(ivorn[, nside]) Creates a MOC object from a given ivorn.
from_json(json_moc) Creates a MOC from a dictionary of HEALPix cell arrays indexed by their depth.
from_lonlat(lon, lat, max_norder) Creates a MOC from astropy lon, lat astropy.units.Quantity.
from_polygon(lon, lat[, inside, max_depth]) Creates a MOC from a polygon
from_polygon_skycoord(skycoord[, inside, …]) Creates a MOC from a polygon.
from_skycoords(skycoords, max_norder) Creates a MOC from an astropy.coordinates.SkyCoord.
from_str(value) Create a MOC from a str.
from_url(url) Creates a MOC object from a given url.
from_vizier_table(table_id[, nside]) Creates a MOC object from a VizieR table.
get_boundaries(self[, order]) Returns the sky coordinates defining the border(s) of the MOC.
intersection(self, another_moc, \*args) Intersection between the MOC instance and other MOCs.
plot(self[, title, frame]) Plot the MOC object using a mollweide projection.
query_simbad(self[, max_rows]) Query a view of SIMBAD data for SIMBAD objects in the coverage of the MOC instance.
query_vizier_table(self, table_id[, max_rows]) Query a VizieR table for sources in the coverage of the MOC instance.
remove_neighbours(self) Removes from the MOC instance the HEALPix cells located at its border.
serialize(self[, format, optional_kw_dict]) Serializes the MOC into a specific format.
union(self, another_moc, \*args) Union between the MOC instance and other MOCs.
write(self, path[, format, overwrite, …]) Writes the MOC to a file.
refine_to_order  
add_neighbours(self)[source]

Extends the MOC instance so that it includes the HEALPix cells touching its border.

The depth of the HEALPix cells added at the border is equal to the maximum depth of the MOC instance.

Returns:
moc : MOC

self extended by one degree of neighbours.

border(self, ax, wcs, **kw_mpl_pathpatch)[source]

Draws the MOC border(s) on a matplotlib axis.

This performs the projection of the sky coordinates defining the perimeter of the MOC to the pixel image coordinate system. You are able to specify various styling kwargs for matplotlib.patches.PathPatch (see the list of valid keywords).

Parameters:
ax : matplotlib.axes.Axes

Matplotlib axis.

wcs : astropy.wcs.WCS

WCS defining the World system <-> Image system projection.

kw_mpl_pathpatch

Plotting arguments for matplotlib.patches.PathPatch

Examples

>>> from mocpy import MOC, WCS
>>> from astropy.coordinates import Angle, SkyCoord
>>> import astropy.units as u
>>> # Load a MOC, e.g. the MOC of GALEXGR6-AIS-FUV
>>> filename = './../resources/P-GALEXGR6-AIS-FUV.fits'
>>> moc = MOC.from_fits(filename)
>>> # Plot the MOC using matplotlib
>>> import matplotlib.pyplot as plt
>>> fig = plt.figure(111, figsize=(15, 15))
>>> # Define a WCS as a context
>>> with WCS(fig, 
...         fov=50 * u.deg,
...         center=SkyCoord(0, 20, unit='deg', frame='icrs'),
...         coordsys="icrs",
...         rotation=Angle(0, u.degree),
...         projection="AIT") as wcs:
...     ax = fig.add_subplot(1, 1, 1, projection=wcs)
...     # Call border giving the matplotlib axe and the `~astropy.wcs.WCS` object.
...     moc.border(ax=ax, wcs=wcs, alpha=0.5, color="red")
>>> plt.xlabel('ra')
>>> plt.ylabel('dec')
>>> plt.grid(color="black", linestyle="dotted")
complement(self)

Returns the complement of the MOC instance.

Returns:
result : MOC or TimeMOC

The resulting MOC.

contains(self, ra, dec, keep_inside=True)[source]

Returns a boolean mask array of the positions lying inside (or outside) the MOC instance.

Parameters:
ra : astropy.units.Quantity

Right ascension array

dec : astropy.units.Quantity

Declination array

keep_inside : bool, optional

True by default. If so the mask describes coordinates lying inside the MOC. If keep_inside is false, contains will return the mask of the coordinates lying outside the MOC.

Returns:
array : ndarray

A mask boolean array

degrade_to_order(self, new_order)

Degrades the MOC instance to a new, less precise, MOC.

The maximum depth (i.e. the depth of the smallest HEALPix cells that can be found in the MOC) of the degraded MOC is set to new_order.

Parameters:
new_order : int

Maximum depth of the output degraded MOC.

Returns:
moc : MOC or TimeMOC

The degraded MOC.

difference(self, another_moc, *args)

Difference between the MOC instance and other MOCs.

Parameters:
another_moc : MOC

The MOC used that will be substracted to self.

args : MOC

Other additional MOCs to perform the difference with.

Returns:
result : MOC or TimeMOC

The resulting MOC.

empty(self)

Checks whether the MOC is empty.

A MOC is empty when its list of HEALPix cell ranges is empty.

Returns:
result: bool

True if the MOC instance is empty.

fill(self, ax, wcs, **kw_mpl_pathpatch)[source]

Draws the MOC on a matplotlib axis.

This performs the projection of the cells from the world coordinate system to the pixel image coordinate system. You are able to specify various styling kwargs for matplotlib.patches.PathPatch (see the list of valid keywords).

Parameters:
ax : matplotlib.axes.Axes

Matplotlib axis.

wcs : astropy.wcs.WCS

WCS defining the World system <-> Image system projection.

kw_mpl_pathpatch

Plotting arguments for matplotlib.patches.PathPatch.

Examples

>>> from mocpy import MOC, WCS
>>> from astropy.coordinates import Angle, SkyCoord
>>> import astropy.units as u
>>> # Load a MOC, e.g. the MOC of GALEXGR6-AIS-FUV
>>> filename = './../resources/P-GALEXGR6-AIS-FUV.fits'
>>> moc = MOC.from_fits(filename)
>>> # Plot the MOC using matplotlib
>>> import matplotlib.pyplot as plt
>>> fig = plt.figure(111, figsize=(15, 15))
>>> # Define a WCS as a context
>>> with WCS(fig, 
...         fov=50 * u.deg,
...         center=SkyCoord(0, 20, unit='deg', frame='icrs'),
...         coordsys="icrs",
...         rotation=Angle(0, u.degree),
...         projection="AIT") as wcs:
...     ax = fig.add_subplot(1, 1, 1, projection=wcs)
...     # Call fill giving the matplotlib axe and the `~astropy.wcs.WCS` object.
...     # We will set the matplotlib keyword linewidth to 0 so that it does not plot
...     # the border of each HEALPix cell.
...     # The color can also be specified along with an alpha value.
...     moc.fill(ax=ax, wcs=wcs, linewidth=0, alpha=0.5, fill=True, color="green")
>>> plt.xlabel('ra')
>>> plt.ylabel('dec')
>>> plt.grid(color="black", linestyle="dotted")
classmethod from_fits(filename)

Loads a MOC from a FITS file.

The specified FITS file must store the MOC (i.e. the list of HEALPix cells it contains) in a binary HDU table.

Parameters:
filename : str

The path to the FITS file.

Returns:
result : MOC or TimeMOC

The resulting MOC.

classmethod from_fits_images(path_l, max_norder)[source]

Loads a MOC from a set of FITS file images.

Parameters:
path_l : [str]

A list of path where the fits image are located.

max_norder : int

The MOC resolution.

Returns:
moc : MOC

The union of all the MOCs created from the paths found in path_l.

classmethod from_healpix_cells(ipix, depth, fully_covered=None)[source]

Creates a MOC from a set of HEALPix cells at a given depth.

Parameters:
ipix : numpy.ndarray

HEALPix cell indices. dtype must be np.uint64

depth : numpy.ndarray

Depth of the HEALPix cells. Must be of the same size of ipix. dtype must be np.uint8

fully_covered : numpy.ndarray, optional

HEALPix cells coverage flags. This flag informs whether a cell is fully covered by a cone (resp. polygon, elliptical cone) or not. Must be of the same size of ipix.

Returns:
moc : MOC

The MOC

Raises:
IndexError

When ipix, depth and fully_covered do not have the same shape

classmethod from_image(header, max_norder, mask=None)[source]

Creates a MOC from an image stored as a FITS file.

Parameters:
header : astropy.io.fits.Header

FITS header containing all the info of where the image is located (position, size, etc…)

max_norder : int

The moc resolution.

mask : numpy.ndarray, optional

A boolean array of the same size of the image where pixels having the value 1 are part of the final MOC and pixels having the value 0 are not.

Returns:
moc : MOC

The resulting MOC.

classmethod from_ivorn(ivorn, nside=256)[source]

Creates a MOC object from a given ivorn.

Parameters:
ivorn : str
nside : int, optional

256 by default

Returns:
result : MOC

The resulting MOC.

classmethod from_json(json_moc)

Creates a MOC from a dictionary of HEALPix cell arrays indexed by their depth.

Parameters:
json_moc : dict(str

A dictionary of HEALPix cell arrays indexed by their depth.

Returns:
moc : MOC or TimeMOC

the MOC.

classmethod from_lonlat(lon, lat, max_norder)[source]

Creates a MOC from astropy lon, lat astropy.units.Quantity.

Parameters:
lon : astropy.units.Quantity

The longitudes of the sky coordinates belonging to the MOC.

lat : astropy.units.Quantity

The latitudes of the sky coordinates belonging to the MOC.

max_norder : int

The depth of the smallest HEALPix cells contained in the MOC.

Returns:
result : MOC

The resulting MOC

classmethod from_polygon(lon, lat, inside=None, max_depth=10)[source]

Creates a MOC from a polygon

The polygon is given as lon and lat astropy.units.Quantity that define the vertices of the polygon. Concave and convex polygons are accepted but self-intersecting ones are currently not properly handled.

Parameters:
lon : astropy.units.Quantity

The longitudes defining the polygon. Can describe convex and concave polygons but not self-intersecting ones.

lat : astropy.units.Quantity

The latitudes defining the polygon. Can describe convex and concave polygons but not self-intersecting ones.

inside : astropy.coordinates.SkyCoord, optional

A point that will be inside the MOC is needed as it is not possible to determine the inside area of a polygon on the unit sphere (there is no infinite area that can be considered as the outside because on the sphere, a closed polygon delimits two finite areas). Possible improvement: take the inside area as the one covering the smallest region on the sphere.

If inside=None (default behavior), the mean of all the vertices is taken as lying inside the polygon. That approach may not work for concave polygons.

max_depth : int, optional

The resolution of the MOC. Set to 10 by default.

Returns:
result : MOC

The resulting MOC

classmethod from_polygon_skycoord(skycoord, inside=None, max_depth=10)[source]

Creates a MOC from a polygon.

The polygon is given as an astropy.coordinates.SkyCoord that contains the vertices of the polygon. Concave and convex polygons are accepted but self-intersecting ones are currently not properly handled.

Parameters:
skycoord : astropy.coordinates.SkyCoord

The sky coordinates defining the vertices of a polygon. It can describe a convex or concave polygon but not a self-intersecting one.

inside : astropy.coordinates.SkyCoord, optional

A point that will be inside the MOC is needed as it is not possible to determine the inside area of a polygon on the unit sphere (there is no infinite area that can be considered as the outside because on the sphere, a closed polygon delimits two finite areas). Possible improvement: take the inside area as the one covering the smallest region on the sphere.

If inside=None (default behavior), the mean of all the vertices is taken as lying inside the polygon. That approach may not work for concave polygons.

max_depth : int, optional

The resolution of the MOC. Set to 10 by default.

Returns:
result : MOC

The resulting MOC

classmethod from_skycoords(skycoords, max_norder)[source]

Creates a MOC from an astropy.coordinates.SkyCoord.

Parameters:
skycoords : astropy.coordinates.SkyCoord

The sky coordinates that will belong to the MOC.

max_norder : int

The depth of the smallest HEALPix cells contained in the MOC.

Returns:
result : MOC

The resulting MOC

classmethod from_str(value)

Create a MOC from a str.

This grammar is expressed is the MOC IVOA specification at section 2.3.2.

Parameters:
value : str

The MOC as a string following the grammar rules.

Returns:
moc : MOC or TimeMOC

The resulting MOC

Examples

>>> from mocpy import MOC
>>> moc = MOC.from_str("2/2-25,28,29 4/0 6/")
classmethod from_url(url)[source]

Creates a MOC object from a given url.

Parameters:
url : str

The url of a FITS file storing a MOC.

Returns:
result : MOC

The resulting MOC.

classmethod from_vizier_table(table_id, nside=256)[source]

Creates a MOC object from a VizieR table.

Info: This method is already implemented in astroquery.cds. You can ask to get a mocpy.moc.MOC object from a vizier catalog ID.

Parameters:
table_id : str

table index

nside : int, optional

256 by default

Returns:
result : MOC

The resulting MOC.

get_boundaries(self, order=None)[source]

Returns the sky coordinates defining the border(s) of the MOC.

The border(s) are expressed as a list of SkyCoord. Each SkyCoord refers to the coordinates of one border of the MOC (i.e. either a border of a connexe MOC part or a border of a hole located in a connexe MOC part). This function is currently not stable: encoding a vertice of a HEALPix cell (N, E, S, W) should not depend on the position of the vertice but rather on the uniq value (+ 2 bits to encode the direction of the vertice).

Parameters:
order : int

The depth of the MOC before computing its boundaries. A shallow depth leads to a faster computation. By default the maximum depth of the MOC is taken.

Raises:
DeprecationWarning

This method is not stable and not tested! A future more stable algorithm will be implemented!

intersection(self, another_moc, *args)

Intersection between the MOC instance and other MOCs.

Parameters:
another_moc : MOC

The MOC used for performing the intersection with self.

args : MOC

Other additional MOCs to perform the intersection with.

Returns:
result : MOC/TimeMOC

The resulting MOC.

max_order

Depth of the smallest HEALPix cells found in the MOC instance.

plot(self, title='MOC', frame=None)[source]

Plot the MOC object using a mollweide projection.

Deprecated: New fill and border methods produce more reliable results and allow you to specify additional matplotlib style parameters.

Parameters:
title : str

The title of the plot

frame : astropy.coordinates.BaseCoordinateFrame, optional

Describes the coordinate system the plot will be (ICRS, Galactic are the only coordinate systems supported).

query_simbad(self, max_rows=10000)[source]

Query a view of SIMBAD data for SIMBAD objects in the coverage of the MOC instance.

query_vizier_table(self, table_id, max_rows=10000)[source]

Query a VizieR table for sources in the coverage of the MOC instance.

remove_neighbours(self)[source]

Removes from the MOC instance the HEALPix cells located at its border.

The depth of the HEALPix cells removed is equal to the maximum depth of the MOC instance.

Returns:
moc : MOC

self minus its HEALPix cells located at its border.

serialize(self, format='fits', optional_kw_dict=None)

Serializes the MOC into a specific format.

Possible formats are FITS, JSON and STRING

Parameters:
format : str

‘fits’ by default. The other possible choice is ‘json’ or ‘str’.

optional_kw_dict : dict

Optional keywords arguments added to the FITS header. Only used if format equals to ‘fits’.

Returns:
result : astropy.io.fits.HDUList or JSON dictionary

The result of the serialization.

sky_fraction

Sky fraction covered by the MOC

union(self, another_moc, *args)

Union between the MOC instance and other MOCs.

Parameters:
another_moc : MOC

The MOC used for performing the union with self.

args : MOC

Other additional MOCs to perform the union with.

Returns:
result : MOC/TimeMOC

The resulting MOC.

write(self, path, format='fits', overwrite=False, optional_kw_dict=None)

Writes the MOC to a file.

Format can be ‘fits’ or ‘json’, though only the fits format is officially supported by the IVOA.

Parameters:
path : str, optional

The path to the file to save the MOC in.

format : str, optional

The format in which the MOC will be serialized before being saved. Possible formats are “fits” or “json”. By default, format is set to “fits”.

overwrite : bool, optional

If the file already exists and you want to overwrite it, then set the overwrite keyword. Default to False.

optional_kw_dict : optional

Optional keywords arguments added to the FITS header. Only used if format equals to ‘fits’.