"""Interpolated frames."""
from __future__ import annotations
__all__ = ["InterpolatedCoordinateFrame", "InterpolatedSkyCoord"]
from typing import TYPE_CHECKING, Any
from astropy.coordinates import (
Angle,
BaseCoordinateFrame,
BaseDifferential,
BaseRepresentation,
BaseRepresentationOrDifferential,
Distance,
SkyCoord,
)
from astropy.utils.decorators import format_doc
from .representation import (
_UNIT_DIF_TYPES,
InterpolatedBaseRepresentationOrDifferential,
InterpolatedRepresentation,
)
from .utils import InterpolatedUnivariateSplinewithUnits as IUSUnits
if TYPE_CHECKING:
from astropy.units import Quantity
from numpy.typing import NDArray
from ._type_hints import CoordinateType, FrameLikeType, RepLikeType
##############################################################################
# CODE
##############################################################################
[docs]
class InterpolatedCoordinateFrame:
"""Wrapper for Coordinate Frame, adding affine interpolations.
Parameters
----------
data : InterpolatedRepresentation or Representation or CoordinateFrame
For either an InterpolatedRepresentation or Representation
the kwarg 'frame' must also be specified.
If CoordinateFrame, then 'frame' is ignored.
affine : Quantity array-like (optional)
if not a Quantity, one is assigned.
Only used if data is not already interpolated.
If data is NOT interpolated, this is required.
Other Parameters
----------------
frame : str or CoordinateFrame
only used if `data` is an InterpolatedRepresentation or Representation
Raises
------
Exception
if `frame` has no error
ValueError
if `data` is not an interpolated type and `affine` is None
TypeError
if `data` is not one of types specified in Parameters.
"""
def __init__(
self,
data: CoordinateType,
affine: Quantity | None = None,
*,
interps: dict[str, Any] | None = None,
**interp_kwargs: Any,
) -> None:
# get rep from CoordinateType
rep = data.data
if isinstance(rep, InterpolatedRepresentation):
pass
elif isinstance(rep, BaseRepresentation):
if affine is None:
msg = "`data` is not already interpolated. Need to pass a Quantity array for `affine`." # noqa: E501
raise ValueError(msg)
rep = InterpolatedRepresentation(rep, affine=affine, interps=interps, **interp_kwargs)
else:
msg = "`data` must be type <InterpolatedRepresentation> or <BaseRepresentation>"
raise TypeError(msg)
self.frame = data.realize_frame(rep)
self._interp_kwargs = interp_kwargs
@property
def _interp_kwargs(self) -> dict[str, Any]:
ikw: dict[str, Any] = self.data._interp_kwargs
return ikw
@_interp_kwargs.setter
def _interp_kwargs(self, value: dict[str, Any]) -> None:
self.data._interp_kwargs = value
[docs]
def __call__(self, affine: Quantity | None = None) -> BaseRepresentation:
"""Evaluate interpolated coordinate frame.
Parameters
----------
affine : `~astropy.units.Quantity` array-like
The affine interpolation parameter.
If None, returns representation points.
Returns
-------
BaseRepresenation
Representation of type ``self.data`` evaluated with `affine`
"""
return self.frame.realize_frame(self.frame.data(affine))
@property
def _class_(self) -> type[InterpolatedCoordinateFrame]:
return object.__class__(self)
def _realize_class(self, data: CoordinateType) -> InterpolatedCoordinateFrame:
return self._class_(data, affine=self.affine, **self._interp_kwargs)
[docs]
def realize_frame(
self,
data: BaseRepresentation,
affine: Quantity | None = None,
**kwargs: Any,
) -> InterpolatedCoordinateFrame:
"""Generates a new frame with new data from another frame.
The new frame may or may not have data. Roughly speaking, the converse
of `replicate_without_data`.
Parameters
----------
data : `~astropy.coordinates.BaseRepresentation`
The representation to use as the data for the new frame.
affine : `~astropy.units.Quantity` | None, optional
The affine interpolation parameter.
**kwargs: Any, optional
Any additional keywords are treated as frame attributes to be set on
the new frame object. In particular, `representation_type` can be
specified.
Returns
-------
frameobj : same as this frame
A new object with the same frame attributes as this one, but with
the ``data`` as the coordinate data.
"""
frame = self.frame.realize_frame(data, **kwargs)
return self._class_(frame, affine=affine, **kwargs)
@property
def uninterpolated(self) -> BaseCoordinateFrame:
"""Return uninterpolated CoordinateFrame."""
return self.frame.realize_frame(self.frame.data.uninterpolated)
#################################################################
# Interpolation Methods
# Mapped to underlying Representation
[docs]
@format_doc(InterpolatedBaseRepresentationOrDifferential.derivative.__doc__) # type: ignore[misc]
def derivative(self, n: int = 1) -> BaseRepresentationOrDifferential:
"""Take nth derivative wrt affine parameter."""
return self.frame.data.derivative(n=n)
@property
def affine(self) -> Quantity: # read-only
return self.frame.data.affine
[docs]
def headless_tangent_vectors(self) -> InterpolatedCoordinateFrame:
r"""Headless tangent vector at each point in affine.
:math:`\vec{x} + \partial_{\affine} \vec{x}(\affine) \Delta\affine`
.. todo::
allow for passing my own points
"""
rep = self.frame.data.headless_tangent_vectors()
return self.realize_frame(rep)
[docs]
def tangent_vectors(self) -> InterpolatedCoordinateFrame:
r"""Tangent vectors along the curve, from the origin.
:math:`\vec{x} + \partial_{\affine} \vec{x}(\affine) \Delta\affine`
.. todo::
allow for passing my own points
"""
rep = self.frame.data.tangent_vectors()
return self.realize_frame(rep)
#################################################################
# Mapping to Underlying CoordinateFrame
@property
def __class__(self) -> type:
"""Make class appear the same as the underlying CoordinateFrame."""
cls: type = self.frame.__class__
return cls
@__class__.setter
def __class__(self, value: Any) -> None: # needed for mypy
msg = "cannot set the `__class__` attribute."
raise TypeError(msg)
def __getattr__(self, key: str) -> Any:
"""Route everything to underlying CoordinateFrame."""
return getattr(object.__getattribute__(self, "frame"), key)
def __len__(self) -> int:
return len(self.frame)
def __getitem__(self, key: int | slice | NDArray) -> InterpolatedCoordinateFrame:
frame = self.frame[key]
affine = self.affine[key]
iframe = self._class_(frame, affine=affine, **self._interp_kwargs)
iframe.representation_type = self.representation_type
return iframe
@property
def representation_type(self) -> BaseRepresentation:
return self.frame.representation_type
@representation_type.setter
def representation_type(self, value: BaseRepresentation) -> None:
self.frame.representation_type = value
[docs]
def represent_as(
self,
base: RepLikeType,
s: str | BaseDifferential = "base",
in_frame_units: bool = False, # noqa: FBT001, FBT002
) -> InterpolatedRepresentation:
"""Generate and return a new representation of this frame's `data`
as a Representation object.
Note: In order to make an in-place change of the representation
of a Frame or SkyCoord object, set the ``representation``
attribute of that object to the desired new representation, or
use the ``set_representation_cls`` method to also set the differential.
Parameters
----------
base : subclass of BaseRepresentation or string
The type of representation to generate. Must be a *class*
(not an instance), or the string name of the representation
class.
s : subclass of `~astropy.coordinates.BaseDifferential`, str, optional
Class in which any velocities should be represented. Must be
a *class* (not an instance), or the string name of the
differential class. If equal to 'base' (default), inferred from
the base class. If `None`, all velocity information is dropped.
in_frame_units : bool, keyword only
Force the representation units to match the specified units
particular to this frame
Returns
-------
newrep : BaseRepresentation-derived object
A new representation object of this frame's `data`.
Raises
------
AttributeError
If this object had no `data`
Examples
--------
>>> from astropy import units as u
>>> from astropy.coordinates import SkyCoord, CartesianRepresentation
>>> coord = SkyCoord(0*u.deg, 0*u.deg)
>>> coord.represent_as(CartesianRepresentation) # doctest: +FLOAT_CMP
<CartesianRepresentation (x, y, z) [dimensionless]
(1., 0., 0.)>
>>> coord.representation_type = CartesianRepresentation
>>> coord # doctest: +FLOAT_CMP
<SkyCoord (ICRS): (x, y, z) [dimensionless]
(1., 0., 0.)>
"""
rep = self.frame.represent_as(base, s=s, in_frame_units=in_frame_units)
return InterpolatedRepresentation(rep, affine=self.affine, **self._interp_kwargs)
[docs]
def copy(self) -> InterpolatedCoordinateFrame:
interp_kwargs = self._interp_kwargs.copy()
frame = self.frame.realize_frame(self.data)
iframe: InterpolatedCoordinateFrame = self._class_(
frame,
affine=self.affine.copy(),
interps=None,
**interp_kwargs,
)
return iframe
def _frame_attrs_repr(self) -> str: # TODO: correct representation
s: str = self.frame._frame_attrs_repr()
return s
def _data_repr(self) -> str: # noqa: C901
"""Returns a string representation of the coordinate data.
This method is modified from the original to include the affine
parameter.
Returns
-------
str
string representation of the data
"""
# if not self.has_data: # must have data to be interpolated
rep_cls = self.representation_type
if rep_cls:
if hasattr(rep_cls, "_unit_representation") and isinstance(
self.frame.data,
rep_cls._unit_representation,
):
rep_cls = self.frame.data.__class__
if "s" in self.frame.data.differentials:
dif_data = self.frame.data.differentials["s"]
dif_cls = (
self.get_representation_cls("s")
if not isinstance(dif_data, _UNIT_DIF_TYPES)
else dif_data.__class__
)
else:
dif_cls = None
data = self.represent_as(rep_cls, dif_cls, in_frame_units=True)
data_repr = repr(data)
# Generate the list of component names out of the repr string
part1, _, remainder = data_repr.partition("(")
if remainder:
comp_str, _, part2 = remainder.partition(")")
comp_names: list[str] = comp_str.split(", ")
affine_name, comp_name_0 = comp_names[0].split("| ")
comp_names[0] = comp_name_0
# Swap in frame-specific component names
rep_comp_names = self.representation_component_names
invnames = {nmrepr: nmpref for nmpref, nmrepr in rep_comp_names.items()}
for i, name in enumerate(comp_names):
comp_names[i] = invnames.get(name, name)
# Reassemble the repr string
data_repr = part1 + "(" + affine_name + "| " + ", ".join(comp_names) + ")" + part2
data_cls_name = "Interpolated" + data.__class__.__name__
if data_repr.startswith("<" + data_cls_name):
# remove both the leading "<" and the space after the name, as well
# as the trailing ">"
i = len(data_cls_name) + 2
data_repr = data_repr[i:-1]
if "s" in self.data.differentials:
data_repr_spl = data_repr.split("\n")
if "has differentials" in data_repr_spl[-1]:
diffrepr = repr(data.differentials["s"]).split("\n")
if diffrepr[0].startswith("<"):
diffrepr[0] = " " + " ".join(diffrepr[0].split(" ")[1:])
for frm_nm, rep_nm in self.get_representation_component_names(
"s",
).items():
diffrepr[0] = diffrepr[0].replace(rep_nm, frm_nm)
if diffrepr[-1].endswith(">"):
diffrepr[-1] = diffrepr[-1][:-1]
data_repr_spl[-1] = "\n".join(diffrepr)
data_repr = "\n".join(data_repr_spl)
return data_repr
def __repr__(self) -> str:
frameattrs = self._frame_attrs_repr()
data_repr = self._data_repr()
if frameattrs:
frameattrs = f" ({frameattrs})"
cls_name = self.__class__.__name__
if data_repr:
return f"<Interpolated{cls_name} Coordinate{frameattrs}: {data_repr}>"
return "TODO!"
# -----------------------------------------------------
# Separation
[docs]
def separation(
self,
point: CoordinateType,
*,
interpolate: bool = True,
affine: Quantity | None = None,
) -> Angle | IUSUnits:
"""Computes on-sky separation between this coordinate and another.
.. note::
If the ``other`` coordinate object is in a different frame, it is
first transformed to the frame of this object. This can lead to
unintuitive behavior if not accounted for. Particularly of note is
that ``self.separation(other)`` and ``other.separation(self)`` may
not give the same answer in this case.
Parameters
----------
point : `~astropy.coordinates.BaseCoordinateFrame`
The coordinate to get the separation to.
interpolate : bool, optional keyword-only
Whether to return the separation as an interpolation, or as points.
affine : `~astropy.units.Quantity` array-like or None, optional
The affine interpolation parameter. If None (default), return
angular width evaluated at all "tick" interpolation points.
Returns
-------
sep : `~astropy.coordinates.Angle` or `~interpolated_coordinates.utils.InterpolatedUnivariateSplinewithUnits`
The on-sky separation between this and the ``other`` coordinate.
Notes
-----
The separation is calculated using the Vincenty formula, which
is stable at all locations, including poles and antipodes [1]_.
.. [1] https://en.wikipedia.org/wiki/Great-circle_distance
""" # noqa: E501
return self._separation(point, angular=True, interpolate=interpolate, affine=affine)
[docs]
def separation_3d(
self,
point: CoordinateType,
*,
interpolate: bool = True,
affine: Quantity | None = None,
) -> Distance | IUSUnits:
"""Compute three dimensional separation between coordinates.
For more on how to use this (and related) functionality, see the
examples in :doc:`astropy:coordinates/matchsep`.
Parameters
----------
point : `~astropy.coordinates.BaseCoordinateFrame`
The coordinate to get the separation to.
interpolate : bool, optional keyword-only
Whether to return the separation as an interpolation, or as points.
affine : `~astropy.units.Quantity` array-like or None, optional
The affine interpolation parameter. If None (default), return
angular width evaluated at all "tick" interpolation points.
Returns
-------
sep : `~astropy.coordinates.Distance` | `~interpolated_coordinates.utils.InterpolatedUnivariateSplinewithUnits`
The real-space distance between these two coordinates.
Raises
------
ValueError
If this or the other coordinate do not have distances.
""" # noqa: E501
return self._separation(point, angular=False, interpolate=interpolate, affine=affine)
def _separation(
self,
point: SkyCoord,
angular: bool, # noqa: FBT001
interpolate: bool, # noqa: FBT001
affine: Quantity | None,
) -> Angle | Distance | IUSUnits:
"""Separation helper function."""
affine = self.affine if affine is None else affine
c = self(affine=affine)
seps = getattr(c, "separation" if angular else "separation_3d")(point)
if not interpolate:
return seps
return IUSUnits(affine, seps) # TODO! if 1 point
#####################################################################
[docs]
class InterpolatedSkyCoord(SkyCoord): # type: ignore[misc]
"""Interpolated SkyCoord."""
def __init__(
self,
*args: Any,
affine: Quantity | None = None,
copy: bool = True,
**kwargs: Any,
) -> None:
keys = tuple(kwargs.keys()) # needed b/c pop changes size
interp_kwargs = {k: kwargs.pop(k) for k in keys if k.startswith("interp_")}
super().__init__(*args, copy=copy, **kwargs)
# change frame to InterpolatedCoordinateFrame
if not isinstance(self.frame, InterpolatedCoordinateFrame):
self._sky_coord_frame = InterpolatedCoordinateFrame(
self.frame,
affine=affine,
**interp_kwargs,
)
[docs]
def __call__(self, affine: Quantity | None = None) -> SkyCoord:
"""Evaluate interpolated representation.
Parameters
----------
affine : `~astropy.units.Quantity` array-like
The affine interpolation parameter.
If None, returns representation points.
Returns
-------
`SkyCoord`
CoordinateFrame of type ``self.frame`` evaluated with `affine`
"""
return SkyCoord(self.frame(affine))
@property
def uninterpolated(self) -> SkyCoord:
"""Return uninterpolated SkyCoord."""
return SkyCoord(self.frame.uninterpolated)
# -----------------------------------------------------
# Separation
[docs]
def separation(
self,
point: CoordinateType,
*,
interpolate: bool = True,
affine: Quantity | None = None,
) -> Angle | IUSUnits:
"""Computes on-sky separation between this coordinate and another.
.. note::
If the ``other`` coordinate object is in a different frame, it is
first transformed to the frame of this object. This can lead to
unintuitive behavior if not accounted for. Particularly of note is
that ``self.separation(other)`` and ``other.separation(self)`` may
not give the same answer in this case.
For more on how to use this (and related) functionality, see the
examples in :doc:`astropy:coordinates/matchsep`.
Parameters
----------
point : `~astropy.coordinates.SkyCoord` or `~astropy.coordinates.BaseCoordinateFrame`
The coordinate to get the separation to.
interpolate : bool, optional keyword-only
Whether to return the separation as an interpolation, or as points.
affine : `~astropy.units.Quantity` array-like or None, optional
The affine interpolation parameter. If None (default), return
angular width evaluated at all "tick" interpolation points.
Returns
-------
sep : ~astropy.coordinates.Angle` | `~interpolated_coordinates.utils.InterpolatedUnivariateSplinewithUnits`
The on-sky separation between this and the ``other`` coordinate.
Notes
-----
The separation is calculated using the Vincenty formula, which
is stable at all locations, including poles and antipodes [1]_.
.. [1] https://en.wikipedia.org/wiki/Great-circle_distance
""" # noqa: E501
return self._separation(point, angular=True, interpolate=interpolate, affine=affine)
[docs]
def separation_3d(
self,
point: CoordinateType,
*,
interpolate: bool = True,
affine: Quantity | None = None,
) -> Distance | IUSUnits:
"""Computes three dimensional separation between this coordinate
and another.
For more on how to use this (and related) functionality, see the
examples in :doc:`astropy:coordinates/matchsep`.
Parameters
----------
point : `~astropy.coordinates.SkyCoord` or `~astropy.coordinates.BaseCoordinateFrame`
The coordinate to get the separation to.
interpolate : bool, optional keyword-only
Whether to return the separation as an interpolation, or as points.
affine : `~astropy.units.Quantity` array-like or None, optional
The affine interpolation parameter. If None (default), return
angular width evaluated at all "tick" interpolation points.
Returns
-------
sep : `~astropy.coordinates.Distance` | `~interpolated_coordinates.utils.InterpolatedUnivariateSplinewithUnits`
The real-space distance between these two coordinates.
Raises
------
ValueError
If this or the other coordinate do not have distances.
""" # noqa: E501
return self._separation(point, angular=False, interpolate=interpolate, affine=affine)
def _separation(
self,
point: SkyCoord,
angular: bool, # noqa: FBT001
interpolate: bool, # noqa: FBT001
affine: Quantity | None,
) -> Angle | Distance | IUSUnits:
"""Separation helper function."""
return InterpolatedCoordinateFrame._separation(
self,
point,
angular=angular,
interpolate=interpolate,
affine=affine,
)
# ---------------------------------------------------------------
[docs]
def match_to_catalog_sky(
self,
catalogcoord: CoordinateType,
nthneighbor: int = 1,
) -> tuple[NDArray, Angle, Quantity]:
"""Finds the nearest on-sky matches of this coordinate in a set of
catalog coordinates.
For more on how to use this (and related) functionality, see the
examples in :doc:`astropy:coordinates/matchsep`.
Parameters
----------
catalogcoord : `~astropy.coordinates.SkyCoord` or `~astropy.coordinates.BaseCoordinateFrame`
The base catalog in which to search for matches. Typically this
will be a coordinate object that is an array (i.e.,
``catalogcoord.isscalar == False``)
nthneighbor : int, optional
Which closest neighbor to search for. Typically ``1`` is
desired here, as that is correct for matching one set of
coordinates to another. The next likely use case is ``2``,
for matching a coordinate catalog against *itself* (``1``
is inappropriate because each point will find itself as the
closest match).
Returns
-------
idx : int array
Indices into ``catalogcoord`` to get the matched points for
each of this object's coordinates. Shape matches this
object.
sep2d : `~astropy.coordinates.Angle`
The on-sky separation between the closest match for each
element in this object in ``catalogcoord``. Shape matches
this object.
dist3d : `~astropy.units.Quantity` ['length']
The 3D distance between the closest match for each element
in this object in ``catalogcoord``. Shape matches this
object. Unless both this and ``catalogcoord`` have associated
distances, this quantity assumes that all sources are at a
distance of 1 (dimensionless).
Notes
-----
This method requires `SciPy <https://www.scipy.org/>`_ to be
installed or it will fail.
See Also
--------
astropy.coordinates.match_coordinates_sky
SkyCoord.match_to_catalog_3d
"""
idx: NDArray # pragma: no cover
sep2d: Angle # pragma: no cover
dist3d: Quantity # pragma: no cover
idx, sep2d, dist3d = super().match_to_catalog_sky(
catalogcoord,
nthneighbor=nthneighbor,
) # pragma: no cover
return idx, sep2d, dist3d # pragma: no cover
[docs]
def match_to_catalog_3d(
self,
catalogcoord: CoordinateType,
nthneighbor: int = 1,
) -> tuple[NDArray, Angle, Quantity]:
"""Finds the nearest 3-dimensional matches of this coordinate to a set
of catalog coordinates.
This finds the 3-dimensional closest neighbor, which is only different
from the on-sky distance if ``distance`` is set in this object or the
``catalogcoord`` object.
For more on how to use this (and related) functionality, see the
examples in :doc:`astropy:coordinates/matchsep`.
Parameters
----------
catalogcoord : `~astropy.coordinates.SkyCoord` or `~astropy.coordinates.BaseCoordinateFrame`
The base catalog in which to search for matches. Typically this
will be a coordinate object that is an array (i.e.,
``catalogcoord.isscalar == False``)
nthneighbor : int, optional
Which closest neighbor to search for. Typically ``1`` is
desired here, as that is correct for matching one set of
coordinates to another. The next likely use case is
``2``, for matching a coordinate catalog against *itself*
(``1`` is inappropriate because each point will find
itself as the closest match).
Returns
-------
idx : int array
Indices into ``catalogcoord`` to get the matched points for
each of this object's coordinates. Shape matches this
object.
sep2d : `~astropy.coordinates.Angle`
The on-sky separation between the closest match for each
element in this object in ``catalogcoord``. Shape matches
this object.
dist3d : `~astropy.units.Quantity` ['length']
The 3D distance between the closest match for each element
in this object in ``catalogcoord``. Shape matches this
object.
Notes
-----
This method requires `SciPy <https://www.scipy.org/>`_ to be
installed or it will fail.
See Also
--------
astropy.coordinates.match_coordinates_3d
SkyCoord.match_to_catalog_sky
"""
idx: NDArray # pragma: no cover
sep2d: Angle # pragma: no cover
dist3d: Quantity # pragma: no cover
idx, sep2d, dist3d = super().match_to_catalog_3d(
catalogcoord,
nthneighbor=nthneighbor,
) # pragma: no cover
return idx, sep2d, dist3d # pragma: no cover
[docs]
def search_around_sky(
self,
searcharoundcoords: CoordinateType,
seplimit: Quantity,
) -> tuple[NDArray, NDArray, Angle, Quantity]:
"""Searches for all coordinates in this object around a supplied set of
points within a given on-sky separation.
This is intended for use on `~astropy.coordinates.SkyCoord` objects
with coordinate arrays, rather than a scalar coordinate. For a scalar
coordinate, it is better to use
`~astropy.coordinates.SkyCoord.separation`.
For more on how to use this (and related) functionality, see the
examples in :doc:`astropy:coordinates/matchsep`.
Parameters
----------
searcharoundcoords : coordinate-like
The coordinates to search around to try to find matching points in
this `SkyCoord`. This should be an object with array coordinates,
not a scalar coordinate object.
seplimit : `~astropy.units.Quantity` ['angle']
The on-sky separation to search within.
Returns
-------
idxsearcharound : int array
Indices into ``searcharoundcoords`` that match the
corresponding elements of ``idxself``. Shape matches
``idxself``.
idxself : int array
Indices into ``self`` that match the
corresponding elements of ``idxsearcharound``. Shape matches
``idxsearcharound``.
sep2d : `~astropy.coordinates.Angle`
The on-sky separation between the coordinates. Shape matches
``idxsearcharound`` and ``idxself``.
dist3d : `~astropy.units.Quantity` ['length']
The 3D distance between the coordinates. Shape matches
``idxsearcharound`` and ``idxself``.
Notes
-----
This method requires `SciPy <https://www.scipy.org/>`_ to be
installed or it will fail.
In the current implementation, the return values are always sorted in
the same order as the ``searcharoundcoords`` (so ``idxsearcharound`` is
in ascending order). This is considered an implementation detail,
though, so it could change in a future release.
See Also
--------
astropy.coordinates.search_around_sky
SkyCoord.search_around_3d
"""
idxsearch: NDArray # pragma: no cover
idxself: NDArray # pragma: no cover
sep2d: Angle # pragma: no cover
dist3d: Quantity # pragma: no cover
idxsearch, idxself, sep2d, dist3d = super().search_around_sky(
searcharoundcoords,
seplimit,
) # pragma: no cover
return idxsearch, idxself, sep2d, dist3d # pragma: no cover
[docs]
def search_around_3d(
self,
searcharoundcoords: CoordinateType,
distlimit: Quantity,
) -> tuple[NDArray, NDArray, Angle, Quantity]:
"""Searches for all coordinates in this object around a point(s).
Searches for all coordinates in this object around a supplied set of
points within a given 3D radius.
This is intended for use on `~astropy.coordinates.SkyCoord` objects
with coordinate arrays, rather than a scalar coordinate. For a scalar
coordinate, it is better to use
`~astropy.coordinates.SkyCoord.separation_3d`.
For more on how to use this (and related) functionality, see the
examples in :doc:`astropy:coordinates/matchsep`.
Parameters
----------
searcharoundcoords : SkyCoord or `~astropy.coordinates.BaseCoordinateFrame`
The coordinates to search around to try to find matching points in
this `SkyCoord`. This should be an object with array coordinates,
not a scalar coordinate object.
distlimit : `~astropy.units.Quantity` ['length']
The physical radius to search within.
Returns
-------
idxsearcharound : int array
Indices into ``searcharoundcoords`` that match the
corresponding elements of ``idxself``. Shape matches
``idxself``.
idxself : int array
Indices into ``self`` that match the
corresponding elements of ``idxsearcharound``. Shape matches
``idxsearcharound``.
sep2d : `~astropy.coordinates.Angle`
The on-sky separation between the coordinates. Shape matches
``idxsearcharound`` and ``idxself``.
dist3d : `~astropy.units.Quantity` ['length']
The 3D distance between the coordinates. Shape matches
``idxsearcharound`` and ``idxself``.
Notes
-----
This method requires `SciPy <https://www.scipy.org/>`_ to be
installed or it will fail.
In the current implementation, the return values are always sorted in
the same order as the ``searcharoundcoords`` (so ``idxsearcharound`` is
in ascending order). This is considered an implementation detail,
though, so it could change in a future release.
See Also
--------
astropy.coordinates.search_around_3d
SkyCoord.search_around_sky
"""
idxsearch: NDArray # pragma: no cover
idxself: NDArray # pragma: no cover
sep2d: Angle # pragma: no cover
dist3d: Quantity # pragma: no cover
idxsearch, idxself, sep2d, dist3d = super().search_around_3d(
searcharoundcoords,
distlimit,
) # pragma: no cover
return idxsearch, idxself, sep2d, dist3d # pragma: no cover