在天文坐标系中处理速度#

使用速度 SkyCoord#

开始获取具有速度的坐标对象的最佳方法是使用 SkyCoord 界面。

实例#

A SkyCoord 要表示一颗具有测量的径向速度但未知自身运动和距离的恒星,可以创建为:

>>> from astropy.coordinates import SkyCoord
>>> import astropy.units as u
>>> sc = SkyCoord(1*u.deg, 2*u.deg, radial_velocity=20*u.km/u.s)
>>> sc  
<SkyCoord (ICRS): (ra, dec) in deg
    (1., 2.)
 (radial_velocity) in km / s
    (20.,)>
>>> sc.radial_velocity  
<Quantity 20.0 km / s>

SkyCoord 以这种方式创建的对象遵循所有相同的变换规则,并且在变换到其他帧时将正确更新其速度。例如,要确定恒星在银河系坐标中的自行,以及在ICRS中测量到的自行:

>>> sc = SkyCoord(1*u.deg, 2*u.deg, pm_ra_cosdec=.2*u.mas/u.yr, pm_dec=.1*u.mas/u.yr)
>>> sc.galactic  
<SkyCoord (Galactic): (l, b) in deg
  ( 99.63785528, -58.70969293)
(pm_l_cosb, pm_b) in mas / yr
  ( 0.22240398,  0.02316181)>

中有关速度支持的有效操作和限制的更多详细信息 astropy.coordinates (特别是 current accuracy limitations ),请参见下面关于较低级别帧对象中的速度支持的更详细讨论。所有这些规则都适用于 SkyCoord 对象,因为它们直接构建在Frame类的速度功能之上,此处详细说明。

使用速度数据创建帧对象#

坐标系类支持存储和转换速度数据(与位置坐标数据一起)。类似于使用 Representation 类抽象出特定的表示,并允许从(例如,笛卡尔到球面)重新表示,速度数据利用 Differential 类来做同样的事情。(有关微分类的详细信息,请参见 表示的微分和导数 )与位置数据一样,微分(速度)组件的名称取决于特定的坐标系。

大多数帧希望速度数据以两个固有运动分量和/或径向速度的形式出现,因为大多数帧的默认差分是 SphericalCosLatDifferential 班级。当有支撑时,所有的自转组件都从 pm_ 并且,默认情况下,纵向构件预计已经包括 cos(latitude) 期限。例如,用于 ICRS 框架是 (pm_ra_cosdecpm_dec

实例#

要使用速度数据以固有运动组件的形式创建帧对象,请执行以下操作:

>>> from astropy.coordinates import ICRS
>>> ICRS(ra=8.67*u.degree, dec=53.09*u.degree,
...      pm_ra_cosdec=4.8*u.mas/u.yr, pm_dec=-15.16*u.mas/u.yr)  
<ICRS Coordinate: (ra, dec) in deg
    (8.67, 53.09)
 (pm_ra_cosdec, pm_dec) in mas / yr
    (4.8, -15.16)>
>>> ICRS(ra=8.67*u.degree, dec=53.09*u.degree,
...      pm_ra_cosdec=4.8*u.mas/u.yr, pm_dec=-15.16*u.mas/u.yr,
...      radial_velocity=23.42*u.km/u.s)  
<ICRS Coordinate: (ra, dec) in deg
    (8.67, 53.09)
 (pm_ra_cosdec, pm_dec, radial_velocity) in (mas / yr, mas / yr, km / s)
    (4.8, -15.16, 23.42)>

对于 Galactic 帧,名称跟踪经度和纬度名称:

>>> from astropy.coordinates import Galactic
>>> Galactic(l=11.23*u.degree, b=58.13*u.degree,
...          pm_l_cosb=21.34*u.mas/u.yr, pm_b=-55.89*u.mas/u.yr)  
<Galactic Coordinate: (l, b) in deg
    (11.23, 58.13)
 (pm_l_cosb, pm_b) in mas / yr
    (21.34, -55.89)>

与位置数据一样,速度数据必须作为 Quantity 物体。

可以更改预期的差分类以控制框架所需的参数名称。默认情况下,固有运动组件应包含 cos(latitude) ,但可以通过指定 SphericalDifferential 类(而不是默认值 SphericalCosLatDifferential ):

>>> from astropy.coordinates import SphericalDifferential
>>> Galactic(l=11.23*u.degree, b=58.13*u.degree,
...          pm_l=21.34*u.mas/u.yr, pm_b=-55.89*u.mas/u.yr,
...          differential_type=SphericalDifferential)  
<Galactic Coordinate: (l, b) in deg
    (11.23, 58.13)
 (pm_l, pm_b) in mas / yr
    (21.34, -55.89)>

这与指定期望的表示类并行工作,只要差分类与表示兼容。例如,要指定笛卡尔坐标系中的所有坐标和速度分量:

>>> from astropy.coordinates import (CartesianRepresentation,
...                                  CartesianDifferential)
>>> Galactic(u=103*u.pc, v=-11*u.pc, w=93.*u.pc,
...          U=31*u.km/u.s, V=-10*u.km/u.s, W=75*u.km/u.s,
...          representation_type=CartesianRepresentation,
...          differential_type=CartesianDifferential)  
<Galactic Coordinate: (u, v, w) in pc
    (103., -11., 93.)
 (U, V, W) in km / s
    (31., -10., 75.)>

请注意 Galactic frame对于笛卡尔位置和速度分量有特殊的标准名称。对于其他帧,这些只是 x,y,zv_x,v_y,v_z ::

>>> ICRS(x=103*u.pc, y=-11*u.pc, z=93.*u.pc,
...      v_x=31*u.km/u.s, v_y=-10*u.km/u.s, v_z=75*u.km/u.s,
...      representation_type=CartesianRepresentation,
...      differential_type=CartesianDifferential)  
<ICRS Coordinate: (x, y, z) in pc
    (103., -11., 93.)
 (v_x, v_y, v_z) in km / s
    (31., -10., 75.)>

对于具有任何表示形式的速度数据的任何帧,也有一些速记,可以更方便地访问常用格式的底层速度数据。对于具有三维速度数据的任何帧对象,可以通过以下方式访问三维笛卡尔速度:

>>> icrs = ICRS(ra=8.67*u.degree, dec=53.09*u.degree,
...             distance=171*u.pc,
...             pm_ra_cosdec=4.8*u.mas/u.yr, pm_dec=-15.16*u.mas/u.yr,
...             radial_velocity=23.42*u.km/u.s)
>>> icrs.velocity 
<CartesianDifferential (d_x, d_y, d_z) in km / s
    ( 23.03160789,  7.44794505,  11.34587732)>

也有检索单曲的速记 Quantity 包含二维固有运动数据并用于检索径向(视线)速度的对象:

>>> icrs.proper_motion 
<Quantity [  4.8 ,-15.16] mas / yr>
>>> icrs.radial_velocity 
<Quantity 23.42 km / s>

向现有帧对象添加速度#

当您有一个现有的Frame对象(或 SkyCoord ),并且想要具有相同位置但增加了速度的对象。从概念上讲,最直接的方法是将差异对象与 realize_frame

实例#

下面的代码片段完成了一个定义明确的情况,在笛卡尔表示法中已知所需的速度。使用将速度添加到现有帧 realize_frame ::

>>> icrs = ICRS(1*u.deg, 2*u.deg, distance=3*u.kpc)
>>> icrs 
<ICRS Coordinate: (ra, dec, distance) in (deg, deg, kpc)
    (1., 2., 3.)>
>>> vel_to_add = CartesianDifferential(4*u.km/u.s, 5*u.km/u.s, 6*u.km/u.s)
>>> newdata = icrs.data.to_cartesian().with_differentials(vel_to_add)
>>> icrs.realize_frame(newdata) 
<ICRS Coordinate: (ra, dec, distance) in (deg, deg, kpc)
    (1., 2., 3.)
 (pm_ra_cosdec, pm_dec, radial_velocity) in (mas / yr, mas / yr, km / s)
    (0.34662023, 0.41161335, 4.29356031)>

即使没有完整的三维坐标(例如,对于距离未知的对象的径向速度观测),也可以使用类似的机构来增加速度。但是,由于缺少显式的单元信息,它需要一种稍微不同的方法来指定差异:

>>> from astropy.coordinates import RadialDifferential
>>> icrs_no_distance = ICRS(1*u.deg, 2*u.deg)
>>> icrs_no_distance
<ICRS Coordinate: (ra, dec) in deg
    (1., 2.)>
>>> rv_to_add = RadialDifferential(500*u.km/u.s)
>>> data_with_rv = icrs_no_distance.data.with_differentials({'s':rv_to_add})
>>> icrs_no_distance.realize_frame(data_with_rv) 
<ICRS Coordinate: (ra, dec) in deg
    (1., 2.)
 (radial_velocity) in km / s
    (500.,)>

我们可以看到,生成的对象与创建对象时指定径向速度时得到的结果相同:

>>> ICRS(1*u.deg, 2*u.deg, radial_velocity=500*u.km/u.s) 
<ICRS Coordinate: (ra, dec) in deg
    (1., 2.)
 (radial_velocity) in km / s
    (500.,)>

用速度变换帧#

将包含速度数据的坐标帧实例转换为不同的帧(可能同时涉及位置和速度变换)的方法与仅变换位置帧实例的方式完全相同。

例子#

要转换包含速度数据的坐标系,请执行以下操作:

>>> from astropy.coordinates import Galactic
>>> icrs = ICRS(ra=8.67*u.degree, dec=53.09*u.degree,
...             pm_ra_cosdec=4.8*u.mas/u.yr, pm_dec=-15.16*u.mas/u.yr)  
>>> icrs.transform_to(Galactic()) 
<Galactic Coordinate: (l, b) in deg
    (120.38084191, -9.69872044)
 (pm_l_cosb, pm_b) in mas / yr
    (3.78957965, -15.44359693)>

然而,速度分量如何变换的细节取决于从起始帧到所需帧所需的特定变换集(即,通过帧变换图的路径)。如果变换链中的所有帧通过 BaseAffineTransform 子类(即矩阵变换或仿射变换),则这些变换可以显式地应用于速度数据。如果不是这样,速度变换是通过对位置变换的有限差分进行数值计算的。有关这两种方法的详细信息,请参见下面的小节。

仿射转换#

包含旋转和/或原点偏移和/或速度偏移的帧变换通过使用 BaseAffineTransform 子类: StaticMatrixTransformDynamicMatrixTransformAffineTransform .

仅矩阵变换(例如,旋转,如 ICRSGalactic )只能在正常运动数据或全空间三维速度上执行。

实例#

要执行仅矩阵变换:

>>> icrs = ICRS(ra=8.67*u.degree, dec=53.09*u.degree,
...             pm_ra_cosdec=4.8*u.mas/u.yr, pm_dec=-15.16*u.mas/u.yr,
...             radial_velocity=23.42*u.km/u.s)
>>> icrs.transform_to(Galactic())  
<Galactic Coordinate: (l, b) in deg
    (120.38084191, -9.69872044)
 (pm_l_cosb, pm_b, radial_velocity) in (mas / yr, mas / yr, km / s)
    (3.78957965, -15.44359693, 23.42)>

相同的旋转矩阵应用于位置向量和速度向量。任何涉及速度偏移的变换都需要所有三维速度分量(通常也需要指定距离),例如, ICRSLSR ::

>>> from astropy.coordinates import LSR
>>> icrs = ICRS(ra=8.67*u.degree, dec=53.09*u.degree,
...             distance=117*u.pc,
...             pm_ra_cosdec=4.8*u.mas/u.yr, pm_dec=-15.16*u.mas/u.yr,
...             radial_velocity=23.42*u.km/u.s)
>>> icrs.transform_to(LSR())  
<LSR Coordinate (v_bary=(11.1, 12.24, 7.25) km / s): (ra, dec, distance) in (deg, deg, pc)
    (8.67, 53.09, 117.)
 (pm_ra_cosdec, pm_dec, radial_velocity) in (mas / yr, mas / yr, km / s)
    (-24.51315607, -2.67935501, 27.07339176)>

有限差分变换#

有些帧变换不能表示为仿射变换。例如,从 AltAz 帧可以包括大气色散校正,这是固有的非线性。此外,有些帧可以更方便地实现为函数,即使它们可以转换为仿射变换。对于这些框架,可以使用有限差分方法来转换速度。注意,这种方法的实现使得用户定义的帧可以以相同的方式使用它(即,通过定义 FunctionTransformWithFiniteDifference 类型)。

这种有限差分方法实际上结合了变换的两个独立(但很重要)元素:

  • 转换 方向 已经存在于起始帧中的速度向量。也就是说,帧变换有时涉及到重新定位坐标系(例如,旋转),新帧中的速度矢量必须考虑到这一点。有限差分法通过沿速度矢量移动起始帧的位置,并计算目标帧中的偏移量来对此进行建模。

  • 由框架运动引起的“诱导”速度 itself 。例如,从以太阳系重心为中心的坐标系向以地球为中心的坐标系的移动包括完全由于地球围绕重心的运动而产生的速度分量。这是通过在略有不同的时间计算目标帧中开始帧的位置,并计算这些位置之间的差异来实现的。请注意,此步骤取决于假设特定的帧属性表示与诱导速度相关的“时间”。按照惯例,这通常是 obstime 帧属性,尽管它是定义有限差分变换函数时可以设置的选项。

例子#

重要的是要认识到有限差分变换具有由有限差分算法和机器精度设定的固有限制。为了说明这个问题,请考虑AltAz到GCRS(即地心)转换。让我们尝试在GCRS框架下计算从地球100 AU处观测到的径向速度,径向速度为10 km/s:

import numpy as np
from matplotlib import pyplot as plt

from astropy import units as u
from astropy.time import Time
from astropy.coordinates import EarthLocation, AltAz, GCRS

time = Time('J2010') + np.linspace(-1,1,1000)*u.min
location = EarthLocation(lon=0*u.deg, lat=45*u.deg)
aa = AltAz(alt=[45]*1000*u.deg, az=90*u.deg, distance=100*u.au,
           radial_velocity=[10]*1000*u.km/u.s,
           location=location, obstime=time)
gcrs = aa.transform_to(GCRS(obstime=time))
plt.plot_date(time.plot_date, gcrs.radial_velocity.to(u.km/u.s))
plt.ylabel('RV [km/s]')

(png, svg, pdf)

../_images/velocities-1.png

这似乎是合理的:径向速度确实应该非常接近10公里/秒,因为帧不涉及速度偏移。

现在让我们考虑100 千帕 作为距离。在这种情况下,我们期望相同:两个帧中的径向速度基本相同:

time = Time('J2010') + np.linspace(-1,1,1000)*u.min
location = EarthLocation(lon=0*u.deg, lat=45*u.deg)
aa = AltAz(alt=[45]*1000*u.deg, az=90*u.deg, distance=100*u.kpc,
           radial_velocity=[10]*1000*u.km/u.s,
           location=location, obstime=time)
gcrs = aa.transform_to(GCRS(obstime=time))
plt.plot_date(time.plot_date, gcrs.radial_velocity.to(u.km/u.s))
plt.ylabel('RV [km/s]')

(png, svg, pdf)

../_images/velocities-2.png

但是这个结果是无稽之谈,取值范围是-1000到1000公里/S,而不是我们预期的~10公里/S。这里的问题的根源是,机器的精度不足以计算距离上的公里数量级上的千分秒级的差值。因此,直接的有限差分方法不适用于使用缺省值的用例。

可以重写发生有限差分的时间步长。例如::

>>> from astropy.coordinates import frame_transform_graph, AltAz, CIRS
>>> trans = frame_transform_graph.get_transform(AltAz, CIRS).transforms[0]
>>> trans.finite_difference_dt = 1 * u.year
>>> gcrs = aa.transform_to(GCRS(obstime=time))  
>>> trans.finite_difference_dt = 1 * u.second  # return to default

在上面的示例中,只有一个转换步骤 AltAzGCRS . 一般来说,两帧之间可能有多个步骤,或者单个步骤可以在内部执行其他转换。可以使用上下文管理器 impose_finite_difference_dt() 要重写的转换图 finite_difference_dt 对于 all 图上的有限差分变换:

>>> from astropy.coordinates import frame_transform_graph
>>> with frame_transform_graph.impose_finite_difference_dt(1 * u.year):
...     gcrs = aa.transform_to(GCRS(obstime=time))  

但要小心这会 not 在上述情况下提供帮助,速度的相关时间范围是秒。(地球相对于特定方向的速度在一年中发生了巨大的变化。)

Astropy的未来版本将对该算法进行改进,以使结果在数值上更稳定,更适合在这些(不寻常)用例中使用。

径向速度修正#

另外,Astropy支持计算重心或日心径向速度修正。虽然将来这可能是使用上述框架的高级便利功能,但当前的实现是独立的,以确保足够的准确性(请参见 径向速度修正 以及 radial_velocity_correction API文档以获取详细信息)。

例子#

这个例子演示了如何计算这种修正,如果在某个特定时间从Keck天文台观测到一个已知的RA和Dec的物体。如果3 m/s左右的精度足够,则可以将计算出的修正值加到任何观测到的径向速度上,以确定最终的日心径向速度:

>>> from astropy.time import Time
>>> from astropy.coordinates import SkyCoord, EarthLocation
>>> # keck = EarthLocation.of_site('Keck')  # the easiest way... but requires internet
>>> keck = EarthLocation.from_geodetic(lat=19.8283*u.deg, lon=-155.4783*u.deg, height=4160*u.m)
>>> sc = SkyCoord(ra=4.88375*u.deg, dec=35.0436389*u.deg)
>>> barycorr = sc.radial_velocity_correction(obstime=Time('2016-6-4'), location=keck)  
>>> barycorr.to(u.km/u.s)  
<Quantity 20.077135 km / s>
>>> heliocorr = sc.radial_velocity_correction('heliocentric', obstime=Time('2016-6-4'), location=keck)  
>>> heliocorr.to(u.km/u.s)  
<Quantity 20.070039 km / s>

请注意,有几种不同的方法来指定校正选项(例如,位置、观测时间等)。见 radial_velocity_correction 有关详细信息,请参阅文档。

精度 radial_velocity_correction#

修正值由 radial_velocity_correction 使用光学近似 \(v = zc\) (见 光谱(多普勒)当量 了解详情)。可以将修正值添加到任何观测到的径向速度上,以提供精确到约3 m/s的校正。如果需要更精确的校正,则必须注意许多细微之处。

首先,你应该总是使用重心修正,因为重心是一个固定点,重力是恒定的。由于日心不满足这些条件,对日心的修正只适用于低精度的工作。因此,为了提高速度,日心改正在 radial_velocity_correction 不包括由于地球表面的电势引起的重力红移等效应。基于这些原因,重心修正 radial_velocity_correction 应始终用于高精度工作。

厘米/秒级径向速度修正所需的其他注意事项见 Wright & Eastman (2014) . 最重要的是,严格地说,重心修正是, 乘法 ,以便将其应用为:

\[v_t=v_m+v_b+\frac{v_b v_m}{c},\]

哪里 \(v_t\) 才是真正的径向速度, \(v_m\) 是测量的径向速度,并且 \(v_b\) 是否由返回的重心校正 radial_velocity_correction 。如果不以这种方式实施重心校正,将导致3m/s量级的误差。

中的重心改正 radial_velocity_correction is consistent with the IDL implementation 莱特和伊斯特曼(2014)论文的水平为10 mm/S,用于无限距离的震源。我们不包括夏皮罗延迟,也不包括那篇论文的公式28中的光行进时间修正。忽略的术语并不重要,除非您需要优于1厘米/S的精度。如果确实需要该精度,请参见 Wright & Eastmann (2014)