1.. _astropy-coordinates-velocities: 2 3Working with Velocities in Astropy Coordinates 4********************************************** 5 6Using Velocities with ``SkyCoord`` 7================================== 8 9The best way to start getting a coordinate object with velocities is to use the 10|SkyCoord| interface. 11 12Examples 13-------- 14 15.. 16 EXAMPLE START 17 Using SkyCoord to Get Coordinate Objects with Velocities 18 19A |SkyCoord| to represent a star with a measured radial velocity but unknown 20proper motion and distance could be created as:: 21 22 >>> from astropy.coordinates import SkyCoord 23 >>> import astropy.units as u 24 >>> sc = SkyCoord(1*u.deg, 2*u.deg, radial_velocity=20*u.km/u.s) 25 >>> sc # doctest: +FLOAT_CMP 26 <SkyCoord (ICRS): (ra, dec) in deg 27 (1., 2.) 28 (radial_velocity) in km / s 29 (20.,)> 30 >>> sc.radial_velocity # doctest: +FLOAT_CMP 31 <Quantity 20.0 km / s> 32 33|SkyCoord| objects created in this manner follow all of the same transformation 34rules and will correctly update their velocities when transformed to other 35frames. For example, to determine proper motions in Galactic coordinates for 36a star with proper motions measured in ICRS:: 37 38 >>> sc = SkyCoord(1*u.deg, 2*u.deg, pm_ra_cosdec=.2*u.mas/u.yr, pm_dec=.1*u.mas/u.yr) 39 >>> sc.galactic # doctest: +FLOAT_CMP 40 <SkyCoord (Galactic): (l, b) in deg 41 ( 99.63785528, -58.70969293) 42 (pm_l_cosb, pm_b) in mas / yr 43 ( 0.22240398, 0.02316181)> 44 45.. 46 EXAMPLE END 47 48For more details on valid operations and limitations of velocity support in 49`astropy.coordinates` (particularly the :ref:`current accuracy limitations 50<astropy-coordinate-finite-difference-velocities>`), see the more detailed 51discussions below of velocity support in the lower-level frame objects. All 52these same rules apply for |SkyCoord| objects, as they are built directly on top 53of the frame classes' velocity functionality detailed here. 54 55.. _astropy-coordinate-custom-frame-with-velocities: 56 57Creating Frame Objects with Velocity Data 58========================================= 59 60The coordinate frame classes support storing and transforming velocity data 61(alongside the positional coordinate data). Similar to the positional data that 62use the ``Representation`` classes to abstract away the particular 63representation and allow re-representing from (e.g., Cartesian to Spherical), 64the velocity data makes use of ``Differential`` classes to do the 65same. (For more information about the differential classes, see 66:ref:`astropy-coordinates-differentials`.) Also like the positional data, the 67names of the differential (velocity) components depend on the particular 68coordinate frame. 69 70Most frames expect velocity data in the form of two proper motion components 71and/or a radial velocity because the default differential for most frames is the 72`~astropy.coordinates.SphericalCosLatDifferential` class. When supported, the 73proper motion components all begin with ``pm_`` and, by default, the 74longitudinal component is expected to already include the ``cos(latitude)`` 75term. For example, the proper motion components for the ``ICRS`` frame are 76(``pm_ra_cosdec``, ``pm_dec``). 77 78Examples 79-------- 80 81.. 82 EXAMPLE START 83 Creating Frame Objects with Proper Motions 84 85To create frame objects with velocity data in the form of proper motion 86components:: 87 88 >>> from astropy.coordinates import ICRS 89 >>> ICRS(ra=8.67*u.degree, dec=53.09*u.degree, 90 ... pm_ra_cosdec=4.8*u.mas/u.yr, pm_dec=-15.16*u.mas/u.yr) # doctest: +FLOAT_CMP 91 <ICRS Coordinate: (ra, dec) in deg 92 (8.67, 53.09) 93 (pm_ra_cosdec, pm_dec) in mas / yr 94 (4.8, -15.16)> 95 >>> ICRS(ra=8.67*u.degree, dec=53.09*u.degree, 96 ... pm_ra_cosdec=4.8*u.mas/u.yr, pm_dec=-15.16*u.mas/u.yr, 97 ... radial_velocity=23.42*u.km/u.s) # doctest: +FLOAT_CMP 98 <ICRS Coordinate: (ra, dec) in deg 99 (8.67, 53.09) 100 (pm_ra_cosdec, pm_dec, radial_velocity) in (mas / yr, mas / yr, km / s) 101 (4.8, -15.16, 23.42)> 102 103For proper motion components in the ``Galactic`` frame, the names track the 104longitude and latitude names:: 105 106 >>> from astropy.coordinates import Galactic 107 >>> Galactic(l=11.23*u.degree, b=58.13*u.degree, 108 ... pm_l_cosb=21.34*u.mas/u.yr, pm_b=-55.89*u.mas/u.yr) # doctest: +FLOAT_CMP 109 <Galactic Coordinate: (l, b) in deg 110 (11.23, 58.13) 111 (pm_l_cosb, pm_b) in mas / yr 112 (21.34, -55.89)> 113 114Like the positional data, velocity data must be passed in as 115`~astropy.units.Quantity` objects. 116 117.. 118 EXAMPLE END 119 120.. 121 EXAMPLE START 122 Changing the Differential Class when Creating Frame Objects 123 124The expected differential class can be changed to control the argument names 125that the frame expects. By default the proper motion components are expected to 126contain the ``cos(latitude)``, but this can be changed by specifying the 127`~astropy.coordinates.SphericalDifferential` class (instead of the default 128`~astropy.coordinates.SphericalCosLatDifferential`):: 129 130 >>> from astropy.coordinates import SphericalDifferential 131 >>> Galactic(l=11.23*u.degree, b=58.13*u.degree, 132 ... pm_l=21.34*u.mas/u.yr, pm_b=-55.89*u.mas/u.yr, 133 ... differential_type=SphericalDifferential) # doctest: +FLOAT_CMP 134 <Galactic Coordinate: (l, b) in deg 135 (11.23, 58.13) 136 (pm_l, pm_b) in mas / yr 137 (21.34, -55.89)> 138 139This works in parallel to specifying the expected representation class, as long 140as the differential class is compatible with the representation. For example, to 141specify all coordinate and velocity components in Cartesian:: 142 143 >>> from astropy.coordinates import (CartesianRepresentation, 144 ... CartesianDifferential) 145 >>> Galactic(u=103*u.pc, v=-11*u.pc, w=93.*u.pc, 146 ... U=31*u.km/u.s, V=-10*u.km/u.s, W=75*u.km/u.s, 147 ... representation_type=CartesianRepresentation, 148 ... differential_type=CartesianDifferential) # doctest: +FLOAT_CMP 149 <Galactic Coordinate: (u, v, w) in pc 150 (103., -11., 93.) 151 (U, V, W) in km / s 152 (31., -10., 75.)> 153 154Note that the ``Galactic`` frame has special, standard names for Cartesian 155position and velocity components. For other frames, these are just ``x,y,z`` and 156``v_x,v_y,v_z``:: 157 158 >>> ICRS(x=103*u.pc, y=-11*u.pc, z=93.*u.pc, 159 ... v_x=31*u.km/u.s, v_y=-10*u.km/u.s, v_z=75*u.km/u.s, 160 ... representation_type=CartesianRepresentation, 161 ... differential_type=CartesianDifferential) # doctest: +FLOAT_CMP 162 <ICRS Coordinate: (x, y, z) in pc 163 (103., -11., 93.) 164 (v_x, v_y, v_z) in km / s 165 (31., -10., 75.)> 166 167.. 168 EXAMPLE END 169 170.. 171 EXAMPLE START 172 Shorthands for Convenient Access to Velocity Data in Frame Objects 173 174For any frame with velocity data with any representation, there are also 175shorthands that provide easier access to the underlying velocity data in 176commonly needed formats. With any frame object with 3D velocity data, the 3D 177Cartesian velocity can be accessed with:: 178 179 >>> icrs = ICRS(ra=8.67*u.degree, dec=53.09*u.degree, 180 ... distance=171*u.pc, 181 ... pm_ra_cosdec=4.8*u.mas/u.yr, pm_dec=-15.16*u.mas/u.yr, 182 ... radial_velocity=23.42*u.km/u.s) 183 >>> icrs.velocity # doctest: +FLOAT_CMP 184 <CartesianDifferential (d_x, d_y, d_z) in km / s 185 ( 23.03160789, 7.44794505, 11.34587732)> 186 187There are also shorthands for retrieving a single `~astropy.units.Quantity` 188object that contains the two-dimensional proper motion data, and for retrieving 189the radial (line-of-sight) velocity:: 190 191 >>> icrs.proper_motion # doctest: +FLOAT_CMP 192 <Quantity [ 4.8 ,-15.16] mas / yr> 193 >>> icrs.radial_velocity # doctest: +FLOAT_CMP 194 <Quantity 23.42 km / s> 195 196.. 197 EXAMPLE END 198 199Adding Velocities to Existing Frame Objects 200=========================================== 201 202Another use case similar to the above comes up when you have an existing frame 203object (or |SkyCoord|) and want an object with the same position but with 204velocities added. The most conceptually direct way to do this is to 205use the differential objects along with 206`~astropy.coordinates.BaseCoordinateFrame.realize_frame`. 207 208Examples 209-------- 210 211.. 212 EXAMPLE START 213 Adding Velocities to Existing Frame Objects 214 215The following snippet accomplishes a well-defined case where the desired 216velocities are known in the Cartesian representation. To add the velocities to 217the existing frame using 218`~astropy.coordinates.BaseCoordinateFrame.realize_frame`:: 219 220 >>> icrs = ICRS(1*u.deg, 2*u.deg, distance=3*u.kpc) 221 >>> icrs # doctest: +FLOAT_CMP 222 <ICRS Coordinate: (ra, dec, distance) in (deg, deg, kpc) 223 (1., 2., 3.)> 224 >>> vel_to_add = CartesianDifferential(4*u.km/u.s, 5*u.km/u.s, 6*u.km/u.s) 225 >>> newdata = icrs.data.to_cartesian().with_differentials(vel_to_add) 226 >>> icrs.realize_frame(newdata) # doctest: +FLOAT_CMP 227 <ICRS Coordinate: (ra, dec, distance) in (deg, deg, kpc) 228 (1., 2., 3.) 229 (pm_ra_cosdec, pm_dec, radial_velocity) in (mas / yr, mas / yr, km / s) 230 (0.34662023, 0.41161335, 4.29356031)> 231 232A similar mechanism can also be used to add velocities even if full 3D coordinates 233are not available (e.g., for a radial velocity observation of an object where 234the distance is unknown). However, it requires a slightly different way of 235specifying the differentials because of the lack of explicit unit information:: 236 237 >>> from astropy.coordinates import RadialDifferential 238 >>> icrs_no_distance = ICRS(1*u.deg, 2*u.deg) 239 >>> icrs_no_distance 240 <ICRS Coordinate: (ra, dec) in deg 241 (1., 2.)> 242 >>> rv_to_add = RadialDifferential(500*u.km/u.s) 243 >>> data_with_rv = icrs_no_distance.data.with_differentials({'s':rv_to_add}) 244 >>> icrs_no_distance.realize_frame(data_with_rv) # doctest: +FLOAT_CMP 245 <ICRS Coordinate: (ra, dec) in deg 246 (1., 2.) 247 (radial_velocity) in km / s 248 (500.,)> 249 250Which we can see yields an object identical to what you get when you specify a 251radial velocity when you create the object:: 252 253 >>> ICRS(1*u.deg, 2*u.deg, radial_velocity=500*u.km/u.s) # doctest: +FLOAT_CMP 254 <ICRS Coordinate: (ra, dec) in deg 255 (1., 2.) 256 (radial_velocity) in km / s 257 (500.,)> 258 259.. 260 EXAMPLE END 261 262.. _astropy-coordinate-transform-with-velocities: 263 264Transforming Frames with Velocities 265=================================== 266 267Transforming coordinate frame instances that contain velocity data to a 268different frame (which may involve both position and velocity transformations) 269is done exactly the same way as transforming position-only frame instances. 270 271Example 272------- 273 274.. 275 EXAMPLE START 276 Transforming Coordinate Frames with Velocities 277 278To transform a coordinate frame that contains velocity data:: 279 280 >>> from astropy.coordinates import Galactic 281 >>> icrs = ICRS(ra=8.67*u.degree, dec=53.09*u.degree, 282 ... pm_ra_cosdec=4.8*u.mas/u.yr, pm_dec=-15.16*u.mas/u.yr) # doctest: +FLOAT_CMP 283 >>> icrs.transform_to(Galactic()) # doctest: +FLOAT_CMP 284 <Galactic Coordinate: (l, b) in deg 285 (120.38084191, -9.69872044) 286 (pm_l_cosb, pm_b) in mas / yr 287 (3.78957965, -15.44359693)> 288 289.. 290 EXAMPLE END 291 292However, the details of how the velocity components are transformed depends on 293the particular set of transforms required to get from the starting frame to the 294desired frame (i.e., the path taken through the frame transform graph). If all 295frames in the chain of transformations are transformed to each other via 296`~astropy.coordinates.BaseAffineTransform` subclasses (i.e., are matrix 297transformations or affine transformations), then the transformations can be 298applied explicitly to the velocity data. If this is not the case, the velocity 299transformation is computed numerically by finite-differencing the positional 300transformation. See the subsections below for more details about these two 301methods. 302 303Affine Transformations 304---------------------- 305 306Frame transformations that involve a rotation and/or an origin shift and/or 307a velocity offset are implemented as affine transformations using the 308`~astropy.coordinates.BaseAffineTransform` subclasses: 309`~astropy.coordinates.StaticMatrixTransform`, 310`~astropy.coordinates.DynamicMatrixTransform`, and 311`~astropy.coordinates.AffineTransform`. 312 313Matrix-only transformations (e.g., rotations such as 314`~astropy.coordinates.ICRS` to `~astropy.coordinates.Galactic`) can be performed 315on proper-motion-only data or full-space, 3D velocities. 316 317Examples 318^^^^^^^^ 319 320.. 321 EXAMPLE START 322 Affine Frame Transformations 323 324To perform a matrix-only transformation:: 325 326 >>> icrs = ICRS(ra=8.67*u.degree, dec=53.09*u.degree, 327 ... pm_ra_cosdec=4.8*u.mas/u.yr, pm_dec=-15.16*u.mas/u.yr, 328 ... radial_velocity=23.42*u.km/u.s) 329 >>> icrs.transform_to(Galactic()) # doctest: +FLOAT_CMP 330 <Galactic Coordinate: (l, b) in deg 331 (120.38084191, -9.69872044) 332 (pm_l_cosb, pm_b, radial_velocity) in (mas / yr, mas / yr, km / s) 333 (3.78957965, -15.44359693, 23.42)> 334 335The same rotation matrix is applied to both the position vector and the velocity 336vector. Any transformation that involves a velocity offset requires all 3D 337velocity components (which typically require specifying a distance as well), 338for example, `~astropy.coordinates.ICRS` to `~astropy.coordinates.LSR`:: 339 340 >>> from astropy.coordinates import LSR 341 >>> icrs = ICRS(ra=8.67*u.degree, dec=53.09*u.degree, 342 ... distance=117*u.pc, 343 ... pm_ra_cosdec=4.8*u.mas/u.yr, pm_dec=-15.16*u.mas/u.yr, 344 ... radial_velocity=23.42*u.km/u.s) 345 >>> icrs.transform_to(LSR()) # doctest: +FLOAT_CMP 346 <LSR Coordinate (v_bary=(11.1, 12.24, 7.25) km / s): (ra, dec, distance) in (deg, deg, pc) 347 (8.67, 53.09, 117.) 348 (pm_ra_cosdec, pm_dec, radial_velocity) in (mas / yr, mas / yr, km / s) 349 (-24.51315607, -2.67935501, 27.07339176)> 350 351.. 352 EXAMPLE END 353 354.. _astropy-coordinate-finite-difference-velocities: 355 356Finite Difference Transformations 357--------------------------------- 358 359Some frame transformations cannot be expressed as affine transformations. 360For example, transformations from the `~astropy.coordinates.AltAz` frame can 361include an atmospheric dispersion correction, which is inherently nonlinear. 362Additionally, some frames are more conveniently implemented as functions, even 363if they can be cast as affine transformations. For these frames, a finite 364difference approach to transforming velocities is available. Note that this 365approach is implemented such that user-defined frames can use it in 366the same manner (i.e., by defining a transformation of the 367`~astropy.coordinates.FunctionTransformWithFiniteDifference` type). 368 369This finite difference approach actually combines two separate (but important) 370elements of the transformation: 371 372 * Transformation of the *direction* of the velocity vector that already exists 373 in the starting frame. That is, a frame transformation sometimes involves 374 reorienting the coordinate frame (e.g., rotation), and the velocity vector 375 in the new frame must account for this. The finite difference approach 376 models this by moving the position of the starting frame along the velocity 377 vector, and computing this offset in the target frame. 378 * The "induced" velocity due to motion of the frame *itself*. For example, 379 shifting from a frame centered at the solar system barycenter to one 380 centered on the Earth includes a velocity component due entirely to the 381 Earth's motion around the barycenter. This is accounted for by computing 382 the location of the starting frame in the target frame at slightly different 383 times, and computing the difference between those. Note that this step 384 depends on assuming that a particular frame attribute represents a "time" 385 of relevance for the induced velocity. By convention this is typically the 386 ``obstime`` frame attribute, although it is an option that can be set when 387 defining a finite difference transformation function. 388 389Example 390^^^^^^^ 391 392.. 393 EXAMPLE START 394 Transforming Velocity Data Between Frames Using a Finite Difference Scheme 395 396It is important to recognize that the finite difference transformations 397have inherent limits set by the finite difference algorithm and machine 398precision. To illustrate this problem, consider the AltAz to GCRS (i.e., 399geocentric) transformation. Let us try to compute the radial velocity in the 400GCRS frame for something observed from the Earth at a distance of 100 AU with a 401radial velocity of 10 km/s: 402 403.. plot:: 404 :context: reset 405 :include-source: 406 407 import numpy as np 408 from matplotlib import pyplot as plt 409 410 from astropy import units as u 411 from astropy.time import Time 412 from astropy.coordinates import EarthLocation, AltAz, GCRS 413 414 time = Time('J2010') + np.linspace(-1,1,1000)*u.min 415 location = EarthLocation(lon=0*u.deg, lat=45*u.deg) 416 aa = AltAz(alt=[45]*1000*u.deg, az=90*u.deg, distance=100*u.au, 417 radial_velocity=[10]*1000*u.km/u.s, 418 location=location, obstime=time) 419 gcrs = aa.transform_to(GCRS(obstime=time)) 420 plt.plot_date(time.plot_date, gcrs.radial_velocity.to(u.km/u.s)) 421 plt.ylabel('RV [km/s]') 422 423This seems plausible: the radial velocity should indeed be very close to 10 km/s 424because the frame does not involve a velocity shift. 425 426Now let us consider 100 *kiloparsecs* as the distance. In this case we expect 427the same: the radial velocity should be essentially the same in both frames: 428 429.. plot:: 430 :context: 431 :include-source: 432 433 time = Time('J2010') + np.linspace(-1,1,1000)*u.min 434 location = EarthLocation(lon=0*u.deg, lat=45*u.deg) 435 aa = AltAz(alt=[45]*1000*u.deg, az=90*u.deg, distance=100*u.kpc, 436 radial_velocity=[10]*1000*u.km/u.s, 437 location=location, obstime=time) 438 gcrs = aa.transform_to(GCRS(obstime=time)) 439 plt.plot_date(time.plot_date, gcrs.radial_velocity.to(u.km/u.s)) 440 plt.ylabel('RV [km/s]') 441 442But this result is nonsense, with values from -1000 to 1000 km/s instead of the 443~10 km/s we expected. The root of the problem here is that the machine 444precision is not sufficient to compute differences on the order of kilometers 445over distances on the order of kiloparsecs. Hence, the straightforward finite 446difference method will not work for this use case with the default values. 447 448.. testsetup:: 449 450 >>> import numpy as np 451 >>> from astropy.coordinates import EarthLocation, AltAz, GCRS 452 >>> from astropy.time import Time 453 >>> time = Time('J2010') + np.linspace(-1,1,1000) * u.min 454 >>> location = EarthLocation(lon=0*u.deg, lat=45*u.deg) 455 >>> aa = AltAz(alt=[45]*1000*u.deg, az=90*u.deg, distance=100*u.kpc, 456 ... radial_velocity=[10]*1000*u.km/u.s, 457 ... location=location, obstime=time) 458 459It is possible to override the timestep over which the finite difference occurs. 460For example:: 461 462 >>> from astropy.coordinates import frame_transform_graph, AltAz, CIRS 463 >>> trans = frame_transform_graph.get_transform(AltAz, CIRS).transforms[0] 464 >>> trans.finite_difference_dt = 1 * u.year 465 >>> gcrs = aa.transform_to(GCRS(obstime=time)) # doctest: +REMOTE_DATA 466 >>> trans.finite_difference_dt = 1 * u.second # return to default 467 468In the above example, there is exactly one transformation step from 469`~astropy.coordinates.AltAz` to `~astropy.coordinates.GCRS`. In general, there 470may be more than one step between two frames, or the single step may perform 471other transformations internally. One can use the context manager 472:func:`~astropy.coordinates.TransformGraph.impose_finite_difference_dt` for the 473transformation graph to override ``finite_difference_dt`` for *all* 474finite-difference transformations on the graph:: 475 476 >>> from astropy.coordinates import frame_transform_graph 477 >>> with frame_transform_graph.impose_finite_difference_dt(1 * u.year): 478 ... gcrs = aa.transform_to(GCRS(obstime=time)) # doctest: +REMOTE_DATA 479 480But beware that this will *not* help in cases like the above, where the relevant 481timescales for the velocities are seconds. (The velocity of the Earth relative 482to a particular direction changes dramatically over the course of one year.) 483 484.. 485 EXAMPLE END 486 487Future versions of Astropy will improve on this algorithm to make the results 488more numerically stable and practical for use in these (not unusual) use cases. 489 490.. _astropy-coordinates-rv-corrs: 491 492Radial Velocity Corrections 493=========================== 494 495Separately from the above, Astropy supports computing barycentric or 496heliocentric radial velocity corrections. While in the future this may 497be a high-level convenience function using the framework described above, the 498current implementation is independent to ensure sufficient accuracy (see 499:ref:`astropy-coordinates-rv-corrs` and the 500`~astropy.coordinates.SkyCoord.radial_velocity_correction` API docs for 501details). 502 503Example 504------- 505 506.. 507 EXAMPLE START 508 Computing Barycentric or Heliocentric Radial Velocity Corrections 509 510This example demonstrates how to compute this correction if observing some 511object at a known RA and Dec from the Keck observatory at a particular time. If 512a precision of around 3 m/s is sufficient, the computed correction can then be 513added to any observed radial velocity to determine the final heliocentric 514radial velocity:: 515 516 >>> from astropy.time import Time 517 >>> from astropy.coordinates import SkyCoord, EarthLocation 518 >>> # keck = EarthLocation.of_site('Keck') # the easiest way... but requires internet 519 >>> keck = EarthLocation.from_geodetic(lat=19.8283*u.deg, lon=-155.4783*u.deg, height=4160*u.m) 520 >>> sc = SkyCoord(ra=4.88375*u.deg, dec=35.0436389*u.deg) 521 >>> barycorr = sc.radial_velocity_correction(obstime=Time('2016-6-4'), location=keck) # doctest: +REMOTE_DATA 522 >>> barycorr.to(u.km/u.s) # doctest: +REMOTE_DATA +FLOAT_CMP 523 <Quantity 20.077135 km / s> 524 >>> heliocorr = sc.radial_velocity_correction('heliocentric', obstime=Time('2016-6-4'), location=keck) # doctest: +REMOTE_DATA 525 >>> heliocorr.to(u.km/u.s) # doctest: +REMOTE_DATA +FLOAT_CMP 526 <Quantity 20.070039 km / s> 527 528Note that there are a few different ways to specify the options for the 529correction (e.g., the location, observation time, etc.). See the 530`~astropy.coordinates.SkyCoord.radial_velocity_correction` docs for more 531information. 532 533.. 534 EXAMPLE END 535 536Precision of `~astropy.coordinates.SkyCoord.radial_velocity_correction` 537------------------------------------------------------------------------ 538 539The correction computed by `~astropy.coordinates.SkyCoord.radial_velocity_correction` 540uses the optical approximation :math:`v = zc` (see :ref:`astropy-units-doppler-equivalencies` 541for details). The correction can be added to any observed radial velocity 542to provide a correction that is accurate to a level of approximately 3 m/s. 543If you need more precise corrections, there are a number of subtleties of 544which you must be aware. 545 546The first is that you should always use a barycentric correction, as the 547barycenter is a fixed point where gravity is constant. Since the heliocenter 548does not satisfy these conditions, corrections to the heliocenter are only 549suitable for low precision work. As a result, and to increase speed, the 550heliocentric correction in 551`~astropy.coordinates.SkyCoord.radial_velocity_correction` does not include 552effects such as the gravitational redshift due to the potential at the Earth's 553surface. For these reasons, the barycentric correction in 554`~astropy.coordinates.SkyCoord.radial_velocity_correction` should always 555be used for high precision work. 556 557Other considerations necessary for radial velocity corrections at the cm/s 558level are outlined in `Wright & Eastman (2014) <https://ui.adsabs.harvard.edu/abs/2014PASP..126..838W>`_. 559Most important is that the barycentric correction is, strictly speaking, 560*multiplicative*, so that you should apply it as: 561 562.. math:: 563 564 v_t = v_m + v_b + \frac{v_b v_m}{c}, 565 566Where :math:`v_t` is the true radial velocity, :math:`v_m` is the measured 567radial velocity and :math:`v_b` is the barycentric correction returned by 568`~astropy.coordinates.SkyCoord.radial_velocity_correction`. Failure to apply 569the barycentric correction in this way leads to errors of order 3 m/s. 570 571The barycentric correction in `~astropy.coordinates.SkyCoord.radial_velocity_correction` is consistent 572with the `IDL implementation <http://astroutils.astronomy.ohio-state.edu/exofast/barycorr.html>`_ of 573the Wright & Eastmann (2014) paper to a level of 10 mm/s for a source at 574infinite distance. We do not include the Shapiro delay nor the light 575travel time correction from equation 28 of that paper. The neglected terms 576are not important unless you require accuracies of better than 1 cm/s. 577If you do require that precision, see `Wright & Eastmann (2014) <https://ui.adsabs.harvard.edu/abs/2014PASP..126..838W>`_. 578