.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "generated/examples/coordinates/plot_sgr-coordinate-frame.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note Click :ref:`here ` to download the full example code .. rst-class:: sphx-glr-example-title .. _sphx_glr_generated_examples_coordinates_plot_sgr-coordinate-frame.py: ========================================================== Create a new coordinate class (for the Sagittarius stream) ========================================================== This document describes in detail how to subclass and define a custom spherical coordinate frame, as discussed in :ref:`astropy:astropy-coordinates-design` and the docstring for `~astropy.coordinates.BaseCoordinateFrame`. In this example, we will define a coordinate system defined by the plane of orbit of the Sagittarius Dwarf Galaxy (hereafter Sgr; as defined in Majewski et al. 2003). The Sgr coordinate system is often referred to in terms of two angular coordinates, :math:`\Lambda,B`. To do this, we need to define a subclass of `~astropy.coordinates.BaseCoordinateFrame` that knows the names and units of the coordinate system angles in each of the supported representations. In this case we support `~astropy.coordinates.SphericalRepresentation` with "Lambda" and "Beta". Then we have to define the transformation from this coordinate system to some other built-in system. Here we will use Galactic coordinates, represented by the `~astropy.coordinates.Galactic` class. See Also -------- * The `gala package `_, which defines a number of Astropy coordinate frames for stellar stream coordinate systems. * Majewski et al. 2003, "A Two Micron All Sky Survey View of the Sagittarius Dwarf Galaxy. I. Morphology of the Sagittarius Core and Tidal Arms", https://arxiv.org/abs/astro-ph/0304198 * Law & Majewski 2010, "The Sagittarius Dwarf Galaxy: A Model for Evolution in a Triaxial Milky Way Halo", https://arxiv.org/abs/1003.1132 * David Law's Sgr info page https://www.stsci.edu/~dlaw/Sgr/ *By: Adrian Price-Whelan, Erik Tollerud* *License: BSD* .. GENERATED FROM PYTHON SOURCE LINES 43-45 Make `print` work the same in all versions of Python, set up numpy, matplotlib, and use a nicer set of plot parameters: .. GENERATED FROM PYTHON SOURCE LINES 45-54 .. code-block:: default import matplotlib.pyplot as plt import numpy as np from astropy.visualization import astropy_mpl_style plt.style.use(astropy_mpl_style) .. GENERATED FROM PYTHON SOURCE LINES 55-56 Import the packages necessary for coordinates .. GENERATED FROM PYTHON SOURCE LINES 56-62 .. code-block:: default import astropy.coordinates as coord import astropy.units as u from astropy.coordinates import frame_transform_graph from astropy.coordinates.matrix_utilities import matrix_transpose, rotation_matrix .. GENERATED FROM PYTHON SOURCE LINES 63-66 The first step is to create a new class, which we'll call ``Sagittarius`` and make it a subclass of `~astropy.coordinates.BaseCoordinateFrame`: .. GENERATED FROM PYTHON SOURCE LINES 66-108 .. code-block:: default class Sagittarius(coord.BaseCoordinateFrame): """ A Heliocentric spherical coordinate system defined by the orbit of the Sagittarius dwarf galaxy, as described in https://ui.adsabs.harvard.edu/abs/2003ApJ...599.1082M and further explained in https://www.stsci.edu/~dlaw/Sgr/. Parameters ---------- representation : `~astropy.coordinates.BaseRepresentation` or None A representation object or None to have no data (or use the other keywords) Lambda : `~astropy.coordinates.Angle`, optional, must be keyword The longitude-like angle corresponding to Sagittarius' orbit. Beta : `~astropy.coordinates.Angle`, optional, must be keyword The latitude-like angle corresponding to Sagittarius' orbit. distance : `~astropy.units.Quantity`, optional, must be keyword The Distance for this object along the line-of-sight. pm_Lambda_cosBeta : `~astropy.units.Quantity`, optional, must be keyword The proper motion along the stream in ``Lambda`` (including the ``cos(Beta)`` factor) for this object (``pm_Beta`` must also be given). pm_Beta : `~astropy.units.Quantity`, optional, must be keyword The proper motion in Declination for this object (``pm_ra_cosdec`` must also be given). radial_velocity : `~astropy.units.Quantity`, optional, keyword-only The radial velocity of this object. """ default_representation = coord.SphericalRepresentation default_differential = coord.SphericalCosLatDifferential frame_specific_representation_info = { coord.SphericalRepresentation: [ coord.RepresentationMapping('lon', 'Lambda'), coord.RepresentationMapping('lat', 'Beta'), coord.RepresentationMapping('distance', 'distance')] } .. GENERATED FROM PYTHON SOURCE LINES 109-127 Breaking this down line-by-line, we define the class as a subclass of `~astropy.coordinates.BaseCoordinateFrame`. Then we include a descriptive docstring. The final lines are class-level attributes that specify the default representation for the data, default differential for the velocity information, and mappings from the attribute names used by representation objects to the names that are to be used by the ``Sagittarius`` frame. In this case we override the names in the spherical representations but don't do anything with other representations like cartesian or cylindrical. Next we have to define the transformation from this coordinate system to some other built-in coordinate system; we will use Galactic coordinates. We can do this by defining functions that return transformation matrices, or by simply defining a function that accepts a coordinate and returns a new coordinate in the new system. Because the transformation to the Sagittarius coordinate system is just a spherical rotation from Galactic coordinates, we'll just define a function that returns this matrix. We'll start by constructing the transformation matrix using pre-determined Euler angles and the ``rotation_matrix`` helper function: .. GENERATED FROM PYTHON SOURCE LINES 127-141 .. code-block:: default SGR_PHI = (180 + 3.75) * u.degree # Euler angles (from Law & Majewski 2010) SGR_THETA = (90 - 13.46) * u.degree SGR_PSI = (180 + 14.111534) * u.degree # Generate the rotation matrix using the x-convention (see Goldstein) SGR_MATRIX = ( np.diag([1.,1.,-1.]) @ rotation_matrix(SGR_PSI, "z") @ rotation_matrix(SGR_THETA, "x") @ rotation_matrix(SGR_PHI, "z") ) .. GENERATED FROM PYTHON SOURCE LINES 142-145 Since we already constructed the transformation (rotation) matrix above, and the inverse of a rotation matrix is just its transpose, the required transformation functions are very simple: .. GENERATED FROM PYTHON SOURCE LINES 145-154 .. code-block:: default @frame_transform_graph.transform(coord.StaticMatrixTransform, coord.Galactic, Sagittarius) def galactic_to_sgr(): """ Compute the transformation matrix from Galactic spherical to heliocentric Sgr coordinates. """ return SGR_MATRIX .. GENERATED FROM PYTHON SOURCE LINES 155-162 The decorator ``@frame_transform_graph.transform(coord.StaticMatrixTransform, coord.Galactic, Sagittarius)`` registers this function on the ``frame_transform_graph`` as a coordinate transformation. Inside the function, we simply return the previously defined rotation matrix. We then register the inverse transformation by using the transpose of the rotation matrix (which is faster to compute than the inverse): .. GENERATED FROM PYTHON SOURCE LINES 162-171 .. code-block:: default @frame_transform_graph.transform(coord.StaticMatrixTransform, Sagittarius, coord.Galactic) def sgr_to_galactic(): """ Compute the transformation matrix from heliocentric Sgr coordinates to spherical Galactic. """ return matrix_transpose(SGR_MATRIX) .. GENERATED FROM PYTHON SOURCE LINES 172-177 Now that we've registered these transformations between ``Sagittarius`` and `~astropy.coordinates.Galactic`, we can transform between *any* coordinate system and ``Sagittarius`` (as long as the other system has a path to transform to `~astropy.coordinates.Galactic`). For example, to transform from ICRS coordinates to ``Sagittarius``, we would do: .. GENERATED FROM PYTHON SOURCE LINES 177-182 .. code-block:: default icrs = coord.SkyCoord(280.161732*u.degree, 11.91934*u.degree, frame='icrs') sgr = icrs.transform_to(Sagittarius) print(sgr) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 183-185 Or, to transform from the ``Sagittarius`` frame to ICRS coordinates (in this case, a line along the ``Sagittarius`` x-y plane): .. GENERATED FROM PYTHON SOURCE LINES 185-191 .. code-block:: default sgr = coord.SkyCoord(Lambda=np.linspace(0, 2*np.pi, 128)*u.radian, Beta=np.zeros(128)*u.radian, frame='sagittarius') icrs = sgr.transform_to(coord.ICRS) print(icrs) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 192-193 As an example, we'll now plot the points in both coordinate systems: .. GENERATED FROM PYTHON SOURCE LINES 193-207 .. code-block:: default fig, axes = plt.subplots(2, 1, figsize=(8, 10), subplot_kw={'projection': 'aitoff'}) axes[0].set_title("Sagittarius") axes[0].plot(sgr.Lambda.wrap_at(180*u.deg).radian, sgr.Beta.radian, linestyle='none', marker='.') axes[1].set_title("ICRS") axes[1].plot(icrs.ra.wrap_at(180*u.deg).radian, icrs.dec.radian, linestyle='none', marker='.') plt.show() .. image-sg:: /generated/examples/coordinates/images/sphx_glr_plot_sgr-coordinate-frame_001.png :alt: Sagittarius, ICRS :srcset: /generated/examples/coordinates/images/sphx_glr_plot_sgr-coordinate-frame_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 208-212 This particular transformation is just a spherical rotation, which is a special case of an Affine transformation with no vector offset. The transformation of velocity components is therefore natively supported as well: .. GENERATED FROM PYTHON SOURCE LINES 212-245 .. code-block:: default sgr = coord.SkyCoord(Lambda=np.linspace(0, 2*np.pi, 128)*u.radian, Beta=np.zeros(128)*u.radian, pm_Lambda_cosBeta=np.random.uniform(-5, 5, 128)*u.mas/u.yr, pm_Beta=np.zeros(128)*u.mas/u.yr, frame='sagittarius') icrs = sgr.transform_to(coord.ICRS) print(icrs) fig, axes = plt.subplots(3, 1, figsize=(8, 10), sharex=True) axes[0].set_title("Sagittarius") axes[0].plot(sgr.Lambda.degree, sgr.pm_Lambda_cosBeta.value, linestyle='none', marker='.') axes[0].set_xlabel(r"$\Lambda$ [deg]") axes[0].set_ylabel( fr"$\mu_\Lambda \, \cos B$ [{sgr.pm_Lambda_cosBeta.unit.to_string('latex_inline')}]") axes[1].set_title("ICRS") axes[1].plot(icrs.ra.degree, icrs.pm_ra_cosdec.value, linestyle='none', marker='.') axes[1].set_ylabel( fr"$\mu_\alpha \, \cos\delta$ [{icrs.pm_ra_cosdec.unit.to_string('latex_inline')}]") axes[2].set_title("ICRS") axes[2].plot(icrs.ra.degree, icrs.pm_dec.value, linestyle='none', marker='.') axes[2].set_xlabel("RA [deg]") axes[2].set_ylabel( fr"$\mu_\delta$ [{icrs.pm_dec.unit.to_string('latex_inline')}]") plt.show() .. image-sg:: /generated/examples/coordinates/images/sphx_glr_plot_sgr-coordinate-frame_002.png :alt: Sagittarius, ICRS, ICRS :srcset: /generated/examples/coordinates/images/sphx_glr_plot_sgr-coordinate-frame_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out Out: .. code-block:: none .. rst-class:: sphx-glr-timing **Total running time of the script:** ( 0 minutes 0.693 seconds) .. _sphx_glr_download_generated_examples_coordinates_plot_sgr-coordinate-frame.py: .. only :: html .. container:: sphx-glr-footer :class: sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_sgr-coordinate-frame.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_sgr-coordinate-frame.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_