.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "generated/examples/coordinates/plot_galactocentric-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_galactocentric-frame.py: ======================================================================== Transforming positions and velocities to and from a Galactocentric frame ======================================================================== This document shows a few examples of how to use and customize the `~astropy.coordinates.Galactocentric` frame to transform Heliocentric sky positions, distance, proper motions, and radial velocities to a Galactocentric, Cartesian frame, and the same in reverse. The main configurable parameters of the `~astropy.coordinates.Galactocentric` frame control the position and velocity of the solar system barycenter within the Galaxy. These are specified by setting the ICRS coordinates of the Galactic center, the distance to the Galactic center (the sun-galactic center line is always assumed to be the x-axis of the Galactocentric frame), and the Cartesian 3-velocity of the sun in the Galactocentric frame. We'll first demonstrate how to customize these values, then show how to set the solar motion instead by inputting the proper motion of Sgr A*. Note that, for brevity, we may refer to the solar system barycenter as just "the sun" in the examples below. *By: Adrian Price-Whelan* *License: BSD* .. GENERATED FROM PYTHON SOURCE LINES 32-34 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 34-43 .. 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 44-45 Import the necessary astropy subpackages .. GENERATED FROM PYTHON SOURCE LINES 45-49 .. code-block:: default import astropy.coordinates as coord import astropy.units as u .. GENERATED FROM PYTHON SOURCE LINES 50-53 Let's first define a barycentric coordinate and velocity in the ICRS frame. We'll use the data for the star HD 39881 from the `Simbad `_ database: .. GENERATED FROM PYTHON SOURCE LINES 53-61 .. code-block:: default c1 = coord.SkyCoord(ra=89.014303*u.degree, dec=13.924912*u.degree, distance=(37.59*u.mas).to(u.pc, u.parallax()), pm_ra_cosdec=372.72*u.mas/u.yr, pm_dec=-483.69*u.mas/u.yr, radial_velocity=0.37*u.km/u.s, frame='icrs') .. GENERATED FROM PYTHON SOURCE LINES 62-66 This is a high proper-motion star; suppose we'd like to transform its position and velocity to a Galactocentric frame to see if it has a large 3D velocity as well. To use the Astropy default solar position and motion parameters, we can simply do: .. GENERATED FROM PYTHON SOURCE LINES 66-69 .. code-block:: default gc1 = c1.transform_to(coord.Galactocentric) .. GENERATED FROM PYTHON SOURCE LINES 70-73 From here, we can access the components of the resulting `~astropy.coordinates.Galactocentric` instance to see the 3D Cartesian velocity components: .. GENERATED FROM PYTHON SOURCE LINES 73-76 .. code-block:: default print(gc1.v_x, gc1.v_y, gc1.v_z) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none 30.254684717897074 km / s 171.29916086104885 km / s 18.19390627095307 km / s .. GENERATED FROM PYTHON SOURCE LINES 77-87 The default parameters for the `~astropy.coordinates.Galactocentric` frame are detailed in the linked documentation, but we can modify the most commonly changes values using the keywords ``galcen_distance``, ``galcen_v_sun``, and ``z_sun`` which set the sun-Galactic center distance, the 3D velocity vector of the sun, and the height of the sun above the Galactic midplane, respectively. The velocity of the sun can be specified as an `~astropy.units.Quantity` object with velocity units and is interpreted as a Cartesian velocity, as in the example below. Note that, as with the positions, the Galactocentric frame is a right-handed system (i.e., the Sun is at negative x values) so ``v_x`` is opposite of the Galactocentric radial velocity: .. GENERATED FROM PYTHON SOURCE LINES 87-94 .. code-block:: default v_sun = [11.1, 244, 7.25] * (u.km / u.s) # [vx, vy, vz] gc_frame = coord.Galactocentric( galcen_distance=8*u.kpc, galcen_v_sun=v_sun, z_sun=0*u.pc) .. GENERATED FROM PYTHON SOURCE LINES 95-96 We can then transform to this frame instead, with our custom parameters: .. GENERATED FROM PYTHON SOURCE LINES 96-100 .. code-block:: default gc2 = c1.transform_to(gc_frame) print(gc2.v_x, gc2.v_y, gc2.v_z) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none 28.427958360720748 km / s 169.69916086104888 km / s 17.70831652451455 km / s .. GENERATED FROM PYTHON SOURCE LINES 101-105 It's sometimes useful to specify the solar motion using the `proper motion of Sgr A* `_ instead of Cartesian velocity components. With an assumed distance, we can convert proper motion components to Cartesian velocity components using `astropy.units`: .. GENERATED FROM PYTHON SOURCE LINES 105-110 .. code-block:: default galcen_distance = 8*u.kpc pm_gal_sgrA = [-6.379, -0.202] * u.mas/u.yr # from Reid & Brunthaler 2004 vy, vz = -(galcen_distance * pm_gal_sgrA).to(u.km/u.s, u.dimensionless_angles()) .. GENERATED FROM PYTHON SOURCE LINES 111-113 We still have to assume a line-of-sight velocity for the Galactic center, which we will again take to be 11 km/s: .. GENERATED FROM PYTHON SOURCE LINES 113-122 .. code-block:: default vx = 11.1 * u.km/u.s v_sun2 = u.Quantity([vx, vy, vz]) # List of Quantity -> a single Quantity gc_frame2 = coord.Galactocentric(galcen_distance=galcen_distance, galcen_v_sun=v_sun2, z_sun=0*u.pc) gc3 = c1.transform_to(gc_frame2) print(gc3.v_x, gc3.v_y, gc3.v_z) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none 28.427958360720748 km / s 167.61484955608267 km / s 18.118916793584443 km / s .. GENERATED FROM PYTHON SOURCE LINES 123-128 The transformations also work in the opposite direction. This can be useful for transforming simulated or theoretical data to observable quantities. As an example, we'll generate 4 theoretical circular orbits at different Galactocentric radii with the same circular velocity, and transform them to Heliocentric coordinates: .. GENERATED FROM PYTHON SOURCE LINES 128-149 .. code-block:: default ring_distances = np.arange(10, 25+1, 5) * u.kpc circ_velocity = 220 * u.km/u.s phi_grid = np.linspace(90, 270, 512) * u.degree # grid of azimuths ring_rep = coord.CylindricalRepresentation( rho=ring_distances[:,np.newaxis], phi=phi_grid[np.newaxis], z=np.zeros_like(ring_distances)[:,np.newaxis]) angular_velocity = (-circ_velocity / ring_distances).to(u.mas/u.yr, u.dimensionless_angles()) ring_dif = coord.CylindricalDifferential( d_rho=np.zeros(phi_grid.shape)[np.newaxis]*u.km/u.s, d_phi=angular_velocity[:,np.newaxis], d_z=np.zeros(phi_grid.shape)[np.newaxis]*u.km/u.s ) ring_rep = ring_rep.with_differentials(ring_dif) gc_rings = coord.SkyCoord(ring_rep, frame=coord.Galactocentric) .. GENERATED FROM PYTHON SOURCE LINES 150-154 First, let's visualize the geometry in Galactocentric coordinates. Here are the positions and velocities of the rings; note that in the velocity plot, the velocities of the 4 rings are identical and thus overlaid under the same curve: .. GENERATED FROM PYTHON SOURCE LINES 154-179 .. code-block:: default fig,axes = plt.subplots(1, 2, figsize=(12,6)) # Positions axes[0].plot(gc_rings.x.T, gc_rings.y.T, marker='None', linewidth=3) axes[0].text(-8., 0, r'$\odot$', fontsize=20) axes[0].set_xlim(-30, 30) axes[0].set_ylim(-30, 30) axes[0].set_xlabel('$x$ [kpc]') axes[0].set_ylabel('$y$ [kpc]') # Velocities axes[1].plot(gc_rings.v_x.T, gc_rings.v_y.T, marker='None', linewidth=3) axes[1].set_xlim(-250, 250) axes[1].set_ylim(-250, 250) axes[1].set_xlabel(f"$v_x$ [{(u.km / u.s).to_string('latex_inline')}]") axes[1].set_ylabel(f"$v_y$ [{(u.km / u.s).to_string('latex_inline')}]") fig.tight_layout() plt.show() .. image-sg:: /generated/examples/coordinates/images/sphx_glr_plot_galactocentric-frame_001.png :alt: plot galactocentric frame :srcset: /generated/examples/coordinates/images/sphx_glr_plot_galactocentric-frame_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 180-182 Now we can transform to Galactic coordinates and visualize the rings in observable coordinates: .. GENERATED FROM PYTHON SOURCE LINES 182-197 .. code-block:: default gal_rings = gc_rings.transform_to(coord.Galactic) fig, ax = plt.subplots(1, 1, figsize=(8, 6)) for i in range(len(ring_distances)): ax.plot(gal_rings[i].l.degree, gal_rings[i].pm_l_cosb.value, label=str(ring_distances[i]), marker='None', linewidth=3) ax.set_xlim(360, 0) ax.set_xlabel('$l$ [deg]') ax.set_ylabel(fr'$\mu_l \, \cos b$ [{(u.mas/u.yr).to_string("latex_inline")}]') ax.legend() plt.show() .. image-sg:: /generated/examples/coordinates/images/sphx_glr_plot_galactocentric-frame_002.png :alt: plot galactocentric frame :srcset: /generated/examples/coordinates/images/sphx_glr_plot_galactocentric-frame_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** ( 0 minutes 0.545 seconds) .. _sphx_glr_download_generated_examples_coordinates_plot_galactocentric-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_galactocentric-frame.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_galactocentric-frame.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_