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()

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

border(ax, wcs, **kw_mpl_pathpatch)

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

complement()

Returns the complement of the MOC instance.

contains(ra, dec[, keep_inside])

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

degrade_to_order(new_order)

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

difference(another_moc, *args)

Difference between the MOC instance and other MOCs.

empty()

Checks whether the MOC is empty.

fill(ax, wcs, **kw_mpl_pathpatch)

Draws the MOC on a matplotlib axis.

from_cone(lon, lat, radius, max_depth[, …])

Creates a MOC from a cone.

from_elliptical_cone(lon, lat, a, b, pa, …)

Creates a MOC from an elliptical cone

from_fits(filename)

Loads a MOC from a FITS file.

from_fits_image(hdu, max_norder[, mask])

Creates a MOC from an image stored as 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_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[, max_depth])

Creates a MOC from a polygon

from_polygon_skycoord(skycoord[, max_depth])

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_valued_healpix_cells(uniq, values[, …])

Creates a MOC from a list of uniq associated with values.

from_vizier_table(table_id[, nside])

Creates a MOC object from a VizieR table.

get_boundaries([order])

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

intersection(another_moc, *args)

Intersection between the MOC instance and other MOCs.

order_to_spatial_resolution(order)

Convert a depth to its equivalent spatial resolution.

plot([title, frame])

Plot the MOC object using a mollweide projection.

query_simbad([max_rows])

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

query_vizier_table(table_id[, max_rows])

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

remove_neighbours()

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

serialize([format, optional_kw_dict])

Serializes the MOC into a specific format.

spatial_resolution_to_order(spatial_resolution)

Convert a spatial resolution to a MOC order.

union(another_moc, *args)

Union between the MOC instance and other MOCs.

write(path[, format, overwrite, …])

Writes the MOC to a file.

refine_to_order

add_neighbours()[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
mocMOC

self extended by one degree of neighbours.

border(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
axmatplotlib.axes.Axes

Matplotlib axis.

wcsastropy.wcs.WCS

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

kw_mpl_pathpatch

Plotting arguments for matplotlib.patches.PathPatch

Examples

>>> from mocpy import MOC, World2ScreenMPL
>>> 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 World2ScreenMPL(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()

Returns the complement of the MOC instance.

Returns
resultMOC or TimeMOC

The resulting MOC.

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

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

Parameters
raastropy.coordinates.Longitude or its supertype astropy.units.Quantity

Right ascension array

decastropy.coordinates.Latitude or its supertype astropy.units.Quantity

Declination array

keep_insidebool, 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
arrayndarray

A mask boolean array

degrade_to_order(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_orderint

Maximum depth of the output degraded MOC.

Returns
mocMOC or TimeMOC

The degraded MOC.

difference(another_moc, *args)

Difference between the MOC instance and other MOCs.

Parameters
another_mocMOC

The MOC used that will be substracted to self.

argsMOC

Other additional MOCs to perform the difference with.

Returns
resultMOC or TimeMOC

The resulting MOC.

empty()

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(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
axmatplotlib.axes.Axes

Matplotlib axis.

wcsastropy.wcs.WCS

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

kw_mpl_pathpatch

Plotting arguments for matplotlib.patches.PathPatch.

Examples

>>> from mocpy import MOC, World2ScreenMPL
>>> 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 World2ScreenMPL(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_cone(lon, lat, radius, max_depth, delta_depth=2)[source]

Creates a MOC from a cone.

The cone is centered around the (lon, lat) position with a radius expressed by radius.

Parameters
lonastropy.coordinates.Longitude or its supertype astropy.units.Quantity

The longitude of the center of the cone.

latastropy.coordinates.Latitude or its supertype astropy.units.Quantity

The latitude of the center of the cone.

radiusastropy.coordinates.Angle

The radius angle of the cone.

max_depthint

Maximum HEALPix cell resolution.

delta_depthint, optional

To control the approximation, you can choose to perform the computations at a deeper depth using the depth_delta parameter. The depth at which the computations will be made will therefore be equal to max_depth + depth_delta.

Returns
resultMOC

The resulting MOC

Examples

>>> from mocpy import MOC
>>> import astropy.units as u
>>> from astropy.coordinates import Angle, Longitude, Latitude
>>> moc = MOC.from_cone(
...  lon=Longitude(0 * u.deg),
...  lat=Latitude(0 * u.deg),
...  radius=Angle(10, u.deg),
...  max_depth=10
... )
classmethod from_elliptical_cone(lon, lat, a, b, pa, max_depth, delta_depth=2)[source]

Creates a MOC from an elliptical cone

The ellipse is centered around the (lon, lat) position. a (resp. b) corresponds to the semi-major axis magnitude (resp. semi-minor axis magnitude). pa is expressed as a Angle and defines the position angle of the elliptical cone.

Parameters
lonastropy.coordinates.Longitude or its supertype astropy.units.Quantity

The longitude of the center of the elliptical cone.

latastropy.coordinates.Latitude or its supertype astropy.units.Quantity

The latitude of the center of the elliptical cone.

aastropy.coordinates.Angle

The semi-major axis angle of the elliptical cone.

bastropy.coordinates.Angle

The semi-minor axis angle of the elliptical cone.

paastropy.coordinates.Angle

The position angle (i.e. the angle between the north and the semi-major axis, east-of-north).

max_depthint

Maximum HEALPix cell resolution.

delta_depthint, optional

To control the approximation, you can choose to perform the computations at a deeper depth using the depth_delta parameter. The depth at which the computations will be made will therefore be equal to depth + depth_delta.

Returns
resultMOC

The resulting MOC

Examples

>>> from mocpy import MOC
>>> import astropy.units as u
>>> from astropy.coordinates import Angle, Longitude, Latitude
>>> moc = MOC.from_elliptical_cone(
...  lon=Longitude(0 * u.deg),
...  lat=Latitude(0 * u.deg),
...  a=Angle(10, u.deg),
...  b=Angle(5, u.deg),
...  pa=Angle(0, u.deg),
...  max_depth=10
... )
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
filenamestr

The path to the FITS file.

Returns
resultMOC or TimeMOC

The resulting MOC.

classmethod from_fits_image(hdu, max_norder, mask=None)[source]

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

Parameters
hduHDU object

HDU containing the data of the image

max_norderint

The moc resolution.

masknumpy.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
mocMOC

The resulting MOC.

classmethod from_fits_images(path_l, max_norder)[source]

Loads a MOC from a set of FITS file images.

Assumes the data of the image is stored in the first HDU of the FITS file. Please call from_fits_image for passing another hdu than the first one.

Parameters
path_l[str]

A list of path where the fits image are located.

max_norderint

The MOC resolution.

Returns
mocMOC

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
ipixnumpy.ndarray

HEALPix cell indices in the NESTED notation. dtype must be np.uint64

depthnumpy.ndarray

Depth of the HEALPix cells. Must be of the same size of ipix. dtype must be np.uint8. Corresponds to the level of an HEALPix cell in astropy.healpix.

fully_coverednumpy.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
mocMOC

The MOC

Raises
IndexError

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

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

Creates a MOC object from a given ivorn.

Parameters
ivornstr
nsideint, optional

256 by default

Returns
resultMOC

The resulting MOC.

classmethod from_json(json_moc)

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

Parameters
json_mocdict(str

A dictionary of HEALPix cell arrays indexed by their depth.

Returns
mocMOC or TimeMOC

the MOC.

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

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

Parameters
lonastropy.coordinates.Longitude or its supertype astropy.units.Quantity

The longitudes of the sky coordinates belonging to the MOC.

latastropy.coordinates.Latitude or its supertype astropy.units.Quantity

The latitudes of the sky coordinates belonging to the MOC.

max_norderint

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

Returns
resultMOC

The resulting MOC

classmethod from_polygon(lon, lat, 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, convex and self-intersecting polygons are accepted.

Parameters
lonastropy.coordinates.Longitude or its supertype astropy.units.Quantity

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

latastropy.coordinates.Latitude or its supertype astropy.units.Quantity

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

max_depthint, optional

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

Returns
resultMOC

The resulting MOC

classmethod from_polygon_skycoord(skycoord, 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, convex and self-intersecting polygons are accepted.

Parameters
skycoordastropy.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.

max_depthint, optional

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

Returns
resultMOC

The resulting MOC

classmethod from_skycoords(skycoords, max_norder)[source]

Creates a MOC from an astropy.coordinates.SkyCoord.

Parameters
skycoordsastropy.coordinates.SkyCoord

The sky coordinates that will belong to the MOC.

max_norderint

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

Returns
resultMOC

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
valuestr

The MOC as a string following the grammar rules.

Returns
mocMOC 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
urlstr

The url of a FITS file storing a MOC.

Returns
resultMOC

The resulting MOC.

classmethod from_valued_healpix_cells(uniq, values, max_depth=None, cumul_from=0.0, cumul_to=1.0)[source]

Creates a MOC from a list of uniq associated with values.

HEALPix cells are first sorted by their values. The MOC contains the cells from which the cumulative value is between cumul_from and cumul_to. Cells being on the fence are recursively splitted and added until the depth of the cells is equal to max_norder.

Parameters
uniqnumpy.ndarray

HEALPix cell indices written in uniq. dtype must be np.uint64

valuesnumpy.ndarray

Probabilities associated with each uniq cells. dtype must be np.float64

max_depthint, optional

The max depth of the MOC. If a depth is given, degrade the MOC to this depth before returning it to the user. Otherwise choose as max_depth the depth corresponding to the smallest HEALPix cell found in uniq.

cumul_fromfloat

Cumulative value from which cells will be added to the MOC

cumul_tofloat

Cumulative value to which cells will be added to the MOC

Returns
resultMOC

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_idstr

table index

nsideint, optional

256 by default

Returns
resultMOC

The resulting MOC.

get_boundaries(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
orderint

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.

Returns
coords: [SkyCoord]

A list of SkyCoord each describing one border.

Raises
DeprecationWarning

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

intersection(another_moc, *args)

Intersection between the MOC instance and other MOCs.

Parameters
another_mocMOC

The MOC used for performing the intersection with self.

argsMOC

Other additional MOCs to perform the intersection with.

Returns
resultMOC/TimeMOC

The resulting MOC.

property max_order

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

static order_to_spatial_resolution(order)[source]

Convert a depth to its equivalent spatial resolution.

Parameters
orderint

Spatial depth.

Returns
spatial_resolutionAngle

Spatial resolution.

plot(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
titlestr

The title of the plot

frameastropy.coordinates.BaseCoordinateFrame, optional

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

query_simbad(max_rows=10000)[source]

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

query_vizier_table(table_id, max_rows=10000)[source]

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

remove_neighbours()[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
mocMOC

self minus its HEALPix cells located at its border.

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

Serializes the MOC into a specific format.

Possible formats are FITS, JSON and STRING

Parameters
formatstr

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

optional_kw_dictdict

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

Returns
resultastropy.io.fits.HDUList or JSON dictionary

The result of the serialization.

property sky_fraction

Sky fraction covered by the MOC

static spatial_resolution_to_order(spatial_resolution)[source]

Convert a spatial resolution to a MOC order.

Parameters
spatial_resolutionAngle

Spatial resolution

Returns
orderint

The order corresponding to the spatial resolution

union(another_moc, *args)

Union between the MOC instance and other MOCs.

Parameters
another_mocMOC

The MOC used for performing the union with self.

argsMOC

Other additional MOCs to perform the union with.

Returns
resultMOC/TimeMOC

The resulting MOC.

write(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
pathstr

The path to the file to save the MOC in.

formatstr, 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”.

overwritebool, optional

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

optional_kw_dictoptional

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