4 期许伟维等 : 天体测量法探测系外行星 423 度法为主要的探测方法 ; 2009 年美国国家航空航天局 National Aeronautics and Space Administration, NASA) 建造的开普

第 57 卷第 4 期天文学报 Vol.57 No.4 2016 年 7 月 ACTA ASTRONOMICA SINICA Jul., 2016 doi: 10.15940/j.cnki.0001-5245.2016.04.004 天体测量法探测系外行星许伟维 1,2,3,4 廖新浩 1,3 周永宏 1,3 许雪晴 1,3 1 中国科学院上海天文台上海 200030) 2 上海科技大学物质学院上海 201210) 3 中国科学院行星科学重点实验室上海 200030) 4 中国科学院大学北京 100049) 摘要在目前已发现的系外行星中, 绝大多数是由视向速度法和凌星法探测得到的, 天体测量法仅发现了 1 颗. gaia 卫星数据即将发布, 天体测量法将逐步成为系外行星探测的重要方法之一. 基于天体测量法给出的恒星位置参数序列, 讨论了在求解行星质量和轨道参数时涉及的动力学条件方程计算问题, 给出了具体微分改正公式, 同时也进行了必要的仿真模拟计算. 建立的方法可以较容易地推广到多行星系统. 关键词系外行星, 天体测量学, 天体力学 : 轨道计算与定轨, 方法 : 数据分析中图分类号 : P135; 文献标识码 : A 1 引言系外行星探测和研究是当今国际天文学研究的热点. 1995 年, Mayor 和 Queloz 采用视向速度法在飞马座 51 附近发现了第 1 颗围绕着主序恒星运动的木星质量量级的系外行星飞马座 51 b [1]. 近年来, 随着天文探测技术的进步, 特别是 Kepler 空间望远镜的成功发射和观测, 新的系外行星不断被发现, 其数量明显增多. 根据系外行星网站 exoplanets.eu) 统计, 截止 2015 年 12 月 1 日, 已确认发现 2004 颗系外行星, 隶属于 1269 个行星系统, 包括 498 个多行星系统 ; 此外, 400 多颗系外行星候选体需进一步认证. 在这些已发现的系外行星中, 绝大多数为气态巨行星和热木星约占 70%), 部分为超级地球 20%) 和地球质量大小的行星 2%), 但尚未确认存在生命的类地行星系统. 系外行星的探测方法分为间接探测和直接成像, 间接探测有视向速度法天体测量法脉冲星计时法凌星法微引力透镜法等 [2], 目前实际观测中采用最多的是视向速度法和凌星法. 天体测量法至今只发现了一颗系外行星 HD 176051 b [3]. 2009 年以前, 由于 HARPS High Accuracy Radial velocity Planet Searcher) 的应用 [4], 视向速 2015-12-08 收到原稿, 2016-03-04 收到修改稿国家自然科学基金项目 11133004) 和中国科学院先导 B 项目 XDBO9000000) 资助 1143880832@qq.com xhliao@shao.ac.cn

4 期许伟维等 : 天体测量法探测系外行星 423 度法为主要的探测方法 ; 2009 年美国国家航空航天局 National Aeronautics and Space Administration, NASA) 建造的开普勒太空望远镜发射升空, 凌星法逐渐成为主流方法, 目前该方法已发现超过 1200 颗系外行星 exoplanets.eu). 2013 年 12 月 19 日, 欧洲空间局 European Space Agency, ESA) 研制的 gaia 空间望远镜升空, 其目的是以前所未有的精度对银河系内数以十亿计的恒星进行观测, 测量他们的位置距离和运动, 其观测效率将比同样由 ESA 发射的依巴谷卫星高出数百万倍 [5]. 因此, 人们对于利用 gaia 开展太阳系与系外行星的探测寄予厚望. 天体测量法是通过测量恒星空间位置变化来反演其附近行星的质量和轨道参数的, 相比于目前通常采用的视向速度法和凌星法, 天体测量法有其独特的优势 : 1) 直接确定行星质量. 视向速度法也是通过测量恒星空间位置变化来反演其附近行星的轨道参数的, 但因其观测方程中仅包含了行星质量 m 与行星轨道倾角的正弦 sin i 的乘积项, 所以只能推算出 m sin i, 不能直接确定行星的质量. 天体测量法对应的观测方程不仅包含了 m cos i 项, 而且包含 m 项, 故可以通过观测资料同时直接推算出行星的质量与轨道倾角 [6]. 2) 寻找长周期行星. 根据开普勒第三定律可知, 行星绕恒星运动的轨道半长径 3 次幂与行星轨道周期的平方成正比 ; 同样质量的行星, 行星轨道半长径越长, 对于恒星绕系统质心运动的扰动就越大, 因此越容易被天体测量法分辨出来. 所以天体测量法对长周期行星探测较为敏感. 3) 检测多星系统是否共面. 在已知的系外行星中, 大多数具有较大的轨道偏心率, 其原因可能是由于受到一颗或多颗大质量天体的引力作用 [7], 该作用力使得行星的轨道面不与恒星的赤道面共面, 若在多星系统中, 该作用力将使得各行星轨道不共面 [8]. 研究行星质量和相互夹角对于多星系统的动力学模型动力学稳定以及轨道共振具有重要意义. 4) 太阳系附近寻找类地行星. 利用空间干涉的天体测量能达到亚微角秒级的精度, 对于造成恒星扰动的行星质量检测可达到 m 地球质量 ) 量级, 满足在太阳系附近寻找类地行星的要求 [9]. 对于由 1 颗恒星和 1 颗行星组成的二体系外行星系统, Wright 和 Howard 于 2009 年建立了基于轨道根数的天体测量法对应的观测方程 [10], 其在赤道坐标系下, 包含了视差自行以及行星质量和轨道根数等参数. 通过分析二体运动的解和参数估计理论, 可以给出这些参数的最优解. 稍加分析即可发现, 基于轨道根数的天体测量方法, 当行星轨道偏心率 e 0 时, 因状态转移矩阵奇异而失效, 特别是对于多行星系统, 基于轨道根数的状态方程和状态转移矩阵计算都非常复杂, 不便于实际应用, 因此需要建立简洁可行的确定系外行星轨道的天体测量方法. 本文针对多体系统动力学特性, 提出采用坐标速度描述, 不仅可以避免 e 0 对应的状态转移矩阵奇异问题, 而且可以显著降低整个计算的复杂性. 本文第 2 节对 gaia 测量数据用于系外行星的确定进行了必要的精度分析 ; 第 3 节较全面地总结了基于轨道根数的天体测量方法涉及的理论基础 ; 第 4 节给出了基于坐标速度的天体测量方法对应的观测方程运动方程和状态转移方程 ; 第 5 节列出了仿真模拟计算结果并作了一些分析 ; 第 6 节推导了行星轨道确定时初始估计满足的约束条件, 由此可以显著提高轨道确定微分改正的收敛效率 ; 第 7 节推导了三体问题 1 颗恒星 2 颗行星 ) 的计算公式 ; 最后 1 节给出本文的结论, 同时也讨论存在的问题及其解决途径.

424 天文学报 57 卷 2 天体测量观测精度分析当恒星附近存在行星时, 由于引力作用, 使得恒星与行星绕着它们的质心以相同的角速度作椭圆运动. 假设恒星和行星围绕系统质心运动轨道半长径分别为 a s 和 a p, 则在地面上观察恒星轨道半长径时, 可产生一个角度 2θ [11] 如图 1). θ m ) 1 p a p G M s D 3 m p 4π 2 M 2 3 s 2 p 3 D, 1) 其中, M s 和 m p 分别为恒星和行星的质量, G 为万有引力常数, D 为恒星与观测者之间的距离, p 为二体运动的轨道周期. 角半径 θ 也称作恒星的天体测量特征角, 表征了行星对恒星的扰动大小, 除行星质量外, 还与恒星的质量轨道周期及观测距离有关. 根据目前已发现的系外行星, 表 1 列出了它们的质量与轨道周期分布范围, 以及所绕恒星的质量和与观测距离的分布范围参见网站 exoplanets.eu), 其中行星质量单位为 M J 木星质量 ), 恒星质量单位为 M 太阳质量 ). 针对上述参数分布范围, 本文选取了 5 种模型, 模型中恒星质量皆取为 1 M, 行星质量轨道周期以及观测距离不等, 计算 5 种模型中恒星的天体测量特征角, 见表 2. 图 1 二体运动的轨道示意图 Fig. 1 The orbital sketch of the two-body motion 表 1 已知系外行星的部分参数范围 Table 1 The partial parameters range of the known exoplanets Planet Mass/M J Period/yr Star Mass/M Distance/pc 10 2 60 10 3 500 0.02 4 1 8500 表 2 几种模型下恒星天体测量特征角 Table 2 The astrometric signature of the star in different models Model Planet mass/m J Star Mass/M Period/yr Distance/pc Signature/µas A 1 1 2 15 100 B 1 1 5 15 186 C 2 1 5 100 55 D 5 1 10 1000 22 E 10 1 20 1000 70

4 期许伟维等 : 天体测量法探测系外行星 425 由于 gaia 对于不同星等的恒星观测精度不同, 所以还需分析 gaia 对于系外行星的可观测性. 表 3 给出了对于太阳质量级的恒星, 当其观测距离在 1000 pc 内时的星等及 gaia 对应的位置测量精度, 由表可知 : 在 1000 pc 内, 恒星的视星等 m V 14.75 mag, gaia 的观测精度高于 10 µas, 其值均小于表 2 中各模型的天体测量特征角, 因此 gaia 可进行有效观测. 表 3 太阳质量级恒星在不同距离的视星等及 gaia 对其的观测精度 Table 3 The apparent magnitude of the sunlike stars at different distances and the accuracies of the gaia Distance/pc 10 100 1000 10000 Apparent magnitude/mag 4.75 9.75 14.75 19.75 m V < 10 mag, 4 µas Accuracies of the gaia [12] m V < 15 mag, 10 µas m V < 20 mag, 200 µas 3 基于轨道根数的观测方程利用天体测量法间接探测系外行星时, 通过观测恒星的周期性位置变化来判断其是否存在. 对于 1 颗孤立的恒星, 可通过 5 个天体测量观测值来描述其在天球上的角位置, 分别是赤道坐标系下的赤经 α 0 赤纬 δ 0 在某一特定的坐标系下某一历元时的值, 现在一般都取在国际天球参考系下 J2000.0 的值 ), 恒星自行相应的正交分量 µ α µ δ, 以及视差 ϖ; 绕着恒星作椭圆运动的行星轨道可由 6 个轨道根数描述, 即轨道半长径 a 倾角 i 升交点经度 Ω 近星点经度 ω 偏心率 e 及平近点角 M; 行星质量 m 可由 p 和 a 通过 Kepler 第三定律确定. 因此, 对一个含有 n 颗行星的系统进行轨道拟合时, 需要 5+7n 个参数. 建立如图 2 所示观测示意图, 恒星的运动轨道处于 O-XY Z 坐标系内, O-xyz 系为观测面所在坐标系, z 轴表示观测者视线方向, x 轴指向北极, y 方向由右手法则确定, i Ω ω 为 3 个方位角. 图 2 二体运动的观测示意图 Fig. 2 The sketch of the observation of the two-body motion

426 天文学报 57 卷恒星在空间位置的变化反映了其椭圆轨道运动在观测面上的投影, 由赤经和赤纬两个位置参数决定, 分别用 ˆηt) ˆξt) 表示, 对应的观测方程为 ˆηt) δ 0 + ϱ δ ϖ + µ δ t t 0 ) + ˆx s, 2) ˆξt) α 0 + ϱ α ϖ + µ α t t 0 ) + ŷ s 等式左边为观测值, 右边各项依次为初始值视差改正自行改正以及天体扰动项, 其中, ϱ α ϱ δ 表示视差改正参数 ; ˆx s ŷ s 表示恒星运动在观测面上两个方向的投影, 单位为角度 ; 为了计算和公式表达的方便, 选取适当的计算单位是必要的, 本文分别采取太阳质量 M 太阳赤道半径 A 作为基本质量和长度单位, 取引力常数 G 1, 则相应的时间单位为 A 3 ) 1/2 [T ], 3) GM 故可通过观测距离 D 将 ˆx s ŷ s 归算到系统长度单位下, 分别用 x s y s 表示 x s d ˆx s A y s d. 4) ŷ s A 考虑恒星附近仅存在 1 颗行星情况, 即 n 1, 由天体力学二体问题分析解可将 x s y s 表示为 x s AXt) + F Y t), 5) y s BXt) + GY t) 其中, Xt) cos Et) e Y t) 1 e 2 sin Et) Et) e sin Et) nt + M, 6) A B F G 称为 Thiele-Innes 常数 [13], 可由轨道根数描述, A a cos Ω cos ω sin Ω sin ω cos i) B a sin Ω cos ω + cos Ω sin ω cos i) F a cos Ω sin ω + sin Ω cos ω cos i) G a sin Ω sin ω cos Ω cos ω cos i). 7) 恒星视差与自行改正可通过一定量天体测量数据计算出来, 故在后面的计算中将视差与自行视为已知量. 至此, 我们可将观测方程改写为 η t) AXt) + F Y t). 8) ξ t) BXt) + GY t)

4 期许伟维等 : 天体测量法探测系外行星 427 这里, η t) 和 ξ t) 包含了自行和视差改正项, 已归算为系统长度单位. 假设初始时刻恒星轨道根数和行星质量的估计值为 a i Ω ω e M 及 m, 为方便起见, 用矢量将其表示为 σ σ 1, σ 2, σ 3, σ 4, σ 5, σ 6, σ 7 ). 对上式在恒星初始运动位置处展开, 有 : 4 [ A η t) Xt)+ F ] 7 [ ] Xt) Y t) Y t) σ i + A+ F σ j + σ i1 i σ i σσ σ j5 j σ j σσ [ ] Xt) Y t) 7 A+ F σ 1 + [AXt)+F Y t)] σ 1 σ σσ P i t) σ i +x t) 1 σσ i1 4 [ B ξ t) Xt)+ G ] 7 [ ], Xt) Y t) Y t) σ i + B+ G σ j + σ i1 i σ i σσ σ j5 j σ j σσ [ ] Xt) Y t) 7 B+ G σ 1 + [BXt)+GY t)] σ 1 σ σσ Q i t) σ i +y t) 1 σσ i1 9) 其中, P i t) Q i t) 依次表示对 7 个参数的微分结果. 如果有 k 次观测, 那么 2k 个观测方程可以写成矩阵形式 a η t 1 ) P 1 t 1 )... P 7 t 1 ) i ξ t 1 ) Q 1 t 1 )... Q 7 t 1 ) Ω...... ω, 10) η t k ) P 1 t k )... P 7 t k ) e ξ t k ) Q 1 t k )... Q 7 t k ) M m 其中, ηt) η t) x t) ξt) ξ t) y t). 11) 上述方程 10) 包含 7 个未知量, 因此问题可解至少需要观测 4 次. 将由估计理论首次解得的 σ 0) 叠加到初始时刻近似轨道上, 进行迭代计算, 当第 l 次迭代给出的 σ l) 小于设定的收敛标准 ε 时, 求解完成. 本文统一设定误差限 ε 10 6. 由上述获得的恒星轨道根数与行星质量的估计值, 可以求出恒星初始时刻的坐标 r s 与速度 rs, 由此根据二体关系式, 可以求出行星此时对应的坐标 r p 与速度 rp, r p M s r s m r p M s r s m. 12) 对于多行星系统, 即 n 2, 对应的运动没有分析解, 故上述方法不再适用.

428 天文学报 57 卷 4 基于坐标速度的观测方程在上一节建立的空间直角坐标系中, 令 r s r p 分别表示恒星和行星的位置矢量, M s 和 m 分别代表恒星与行星的质量. 恒星相对于质心的运动方程为 r s Gm 3 r s m + M s ) 2 r s 3 F. 13) 观测方程与轨道根数法对应的恒星观测方程相同, 即 2) 式. 假定 t t 0 时刻恒星位置与速度和行星质量的估计值分别为 r s t 0 ) r s t 0 ) m p, 其近似值即初始估计 ) 为 r r m, 两者相差为 r r m, 将观测方程在 r r m 处展开有 ηt) x s r + x s r r r + x s m m ξt) y s r + y s r r r + y. 14) s m m 令 rs v, 恒星运动和相应的状态转移矩阵可以通过联立求解下述方程组获得, r s v v F ) ) d rs v dt r r ) d v ) ) F rs dt r r s r ) ) d rs v dt v v ) d v ) ) F rs dt v r s v ) ) d rs v dt m m ) d v ) F + ) ) F rs dt m m r s m, 15) 当 t t 0 时, 状态转移方程涉及的初始条件为 r s r s H, 0, r v 这里, H 为 3 阶单位矩阵. v 0, r 5 轨道特征分析与仿真计算 v v H, r s 0 m v 0 m. 16) 目前, gaia 的巡天数据仍未释放, 因此为了验证上述方法的有效性, 只能进行仿真模

4 期许伟维等 : 天体测量法探测系外行星 429 拟计算. 假定 t 0 时刻 1 组恒星的轨道根数和行星质量 σt 0 ), 由该组参数可模拟出 1 组恒星位置的观测数据, 然后用其反解 σt 0 ), 由此可以验证和分析评估轨道确定的方法和精度. 5.1 轨道特征的影响 5.1.1 小偏心率对于小偏心率轨道, 轨道根数法不能进行轨道确定, 但坐标速度法仍有效, 其原因是 e 0 时近星点经度不能确定所致. 具体体现在当 e 0 时, x s y s 可表示为 : x s a [cos Ω cos ω + M) sin Ω sinω + M) cos i], 17) y s a [sin Ω cos ω + M) + cos Ω sinω + M) cos i] 由上式可以看出 x s ω x s M y s ω y s M. 18) 即状态转移矩阵中有关 ω 和 M 的两列相关, 导致状态转移矩阵奇异, 故对于小偏心率的情况, 轨道根数法失效. 5.1.2 轨道倾角对于倾角 i 0 或 180 的轨道, 两种方法都不适用. 轨道根数法中, 由于 x s y s 中包含的 i 项都以 cos i 出现, 故当 i 0 或 180 时, x s y s 对 i 的偏导数都为 0, 具体体现在 A i B i F i G i 0. 19) 坐标速度法中, 当 i 0 时, r s 中不显含 z 分量, 即 z 和 ż 分量. 故有 r s z r s ż 0. 20) 两种情况都导致线性方程组的系数矩阵奇异, 无唯一解. 此外, 天体测量法描述的是恒星 X Y 方向的运动在观测面上的投影, 对于 Z 方向没有直接约束, 这就造成了当 i π i 时, σt 0 ) 一样能满足观测方程. 但此时除倾角项外, 升交点的经度 Ω 和近地点的经度 ω 也存在相应的变化, 这一点我们可从 Thiele-Innes 常数来理解. 我们将其对应关系整理如下 : i, Ω, ω) i, Ω, ω) i, Ω, ω) ±i, Ω + π, ω + π). 21) i, Ω, ω) ±i, Ω π, ω π) 5.2 观测误差的影响观测即有误差, 因此不同的仿真模型应增加相应的观测误差. 以第 2 节中模型 A 为例, m V 5.63 mag < 10 mag, 此时 gaia 的观测精度 θ 4 µas, 若将其视为天体测量观

430 天文学报 57 卷测误差, 对于观测距离 D 15 pc, 该误差 s 的大小为 / ) 4 180 s θ D 4 µas 15 pc 3600 10 6 15 3.086 10 13 )km 8975.879 km. π 将其归算到系统长度单位 A, 即得该模型的观测误差 s 为 s s 8975.879 km A 696300 km 0.013A. 通过 Matlab 软件可计算出任意两组满足均值为 0 方差为 0.013 的高斯白噪声, 将其分别增加到恒星的 x 和 y 方向的观测数据上, 则可获得含有观测误差的模拟观测数据. 由含有误差观测数据可确定初始根数 σt 0 ) 的估计值 σt 0 ), 两者之差为 σ. 表 4 列出了几种模型对于不同观测误差和不同观测次数时确定的初始轨道估计值误差情况. 模型 A 表示一个轨道周期为 2 yr 的木星绕一个距离地球 15 pc 的太阳作椭圆运动见表 2), 此太阳在天球上位置变化最大角度特征角 ) 为 100 µas, 考虑不同观测误差与观测次数的组合, 形成了 3 组模拟观测数据 : A-1 A-2 和 A-3. 模型 A 的初始轨道根数取为 : a 1.587 au, e 0.4000, i 40.00, Ω 80.00, ω 50.00, M 50.00, m 9.552 10 4 M s ; 轨道确定误差统计分为 x y 方向的残差及轨道根数的相对误差. 以 A 模型为例, 对于特征角为 100 µas 的系外行星系统, 考虑观测随机误差分别为 2 µas 和 4 µas 以及观测点数分别为 20 和 40 的情况. 模型 A-1 对应观测误差 2 µas, 观测点数 20, 由坐标速度法确定的轨道残差的 x 和 y 分量分别为 23.75 µas 和 24.20 µas, 由轨道根数法确定的轨道半长径和轨道偏心率的相对误差分别为 1.65% 和 0.92%, 其他根数的误差均小于 10%, 质量确定的精度较高, 可好于千分之一 ; 模型 A-2 对应观测误差 2 µas, 观测点数 40, 相较于 A-1 各根数的相对误差均显著减小, 即拟合结果更优 ; 模型 A-3 对应观测误差 4 µas, 观测点数 40, 相较于 A-2 大部分根数的相对误差有所增大. 通过比较不同模型对应的计算结果表明, 观测数据越多, 相对误差越小, 轨道确定的拟合结果越可靠. 表 4 模型 A-1 A-2 和 A-3 的计算结果 Table 4 The calculations of the models A-1, A-2, and A-3 Residuals of x and y components and relative Model θ Error /µas /µas Number errors of orbit elements x/µas y/µas a/a e i/180 Ω + ω + M) m/%) /%) /%) /%) /360 /%) A-1 100 2 20 23.75 24.20 1.65 0.92 6.82 3.78 3.06 10 2 A-2 100 2 40 18.79 17.84 0.47 0.24 2.26 0.36 7.20 10 3 A-3 100 4 40 26.51 26.66 0.18 2.67 2.62 0.27 9.40 10 3 5.3 初值的选取由于观测方程 2) 中关于恒星初始轨道根数或恒星初始坐标速度函数的表示均是非线性的, 故在轨道确定中需要采用微分迭代方法求解轨道估计值, 因此需要选定恰当的迭代初值, 以使求解过程收敛. 对于轨道初值的选取, 最简单的方法是在一定区间内对

4 期许伟维等 : 天体测量法探测系外行星 431 每个待解参数进行扫描计算, 当迭代初值被选取在参数收敛区间时, 便可给出估计值. 本文讨论的问题涉及的参数收敛区间是 7 维的, 因此要获得收敛的估计值, 需要较大的扫描计算量. 为了方便计算, 采用固定 Ω ω M m 进行 a i e 3 维检索, 再固定 a i e m 进行 Ω ω M 3 维检索. 计算结果列于表 5, 表中列出的是相对误差百分比值, 需要说明的是其中角度量的误差相对量取 180, 偏心率的误差相对量取为 1, 主要是为了避免小偏心率作为相对量导致相对误差失真. 模型 A-1 和 A-3 的轨道半长径收敛区间约为 5%, 其他根数的收敛区间相对较大, 一般大于 10%, 有些甚至可达 20%. 从初步的扫描计算结果看, 轨道半长径 a 的迭代初值选取是影响计算收敛与否的关键因素, 需要尽量靠近估计值. Model 表 5 模型 A-1 和 A-3 的收敛区间 Table 5 The convergence interval of the models A-1 and A-3 Convergence interval of each orbit element/%) a/a i/180 Ω/180 ω/180 e/1 M/180 A-1 3.0 3.0 2.0 11.0 14.0 14.0 16.0 16.0 15.0 15.0 11.0 11.0 A-3 6.0 6.0 9.0 9.0 6.0 11.0 22.0 22.0 12.0 12.0 13.0 22.0 6 关于初始估计取值的约束在上节讨论中, 为减少计算时间, 对 7 维参数空间仅进行了两个 3 维参数空间扫描, 且行星质量均是固定的, 没有参与扫描计算. 事实上, 可通过恒星空间位置观测资料时间序列对上述 7 维参数空间进行适当的约束, 降低参数空间的维数. 具体做法如下 : 已知由天体测量方法给出的某颗恒星的位置时间序列 x s t j ), y s t j ), j 1, 2,, k, 其包含了运动周期信息, 通过对其傅里叶谱分析 FFT) 可获得轨道周期 p, 再根据开普勒第三定律可以给出恒星轨道半长径 a 和行星质量 m p 之间的关系式, p 2 4π2 µ r3 µ G M s + m p ) m p a r M s + m p, 22) 一般认为 M s 的值可由天体物理方法获得. 傅里叶谱分析通常针对等间距数据均匀采样 ) 进行, 故需要先对恒星的观测数据进行插值, 再利用插值数据进行谱分析. 此外, 还采用了另一种方法 LSQ 最小二乘拟合 ), 将插值后的数据进行最小二乘, 拟合成一条正弦曲线, 该正弦曲线的周期即为恒星运动周期. 已知 gaia 的观测年限为 5 yr, 故进行傅里叶谱分析所选取的模型周期不应该大于 5 yr, 否则无法观测到一个完整周期的数据 ; 又根据奈奎斯特取样定理, 所选模型周期应限制在 2.5 yr 内. 此处, 我们选取 5.2 节中所建立的模型进行计算, 图 3 给出了 A-3 模型任意选取 )x y 方向误差数据的分布及进行傅里叶分析的结果, 并将计算结果列于

432 天文学报 57 卷表 6 中. 表 6 模型 A-1 A-2 A-3 进行最小二乘拟合和傅里叶分析的结果 Table 6 The LSQ and FFT results of the models A-1, A-2, and A-3 Model LSQ/d FFT/d Periods/d A-1 698.49 720.00 730.00 735.77 A-2 712.23 716.67 730.00 744.24 A-3 692.58 716.67 730.00 755.88 图 3 模型 A-3 的 x y 方向数据进行傅里叶分析的结果 Fig. 3 The results of x and y components of the model A-3 with FFT 根据统计, LSQ 和 FFT 均能实现观测数据的周期分析. 模型 A 的精确周期为 730 d, 3 个子模型进行 FFT 的结果分别为 720 d 716.67 d 和 716.67 d, 明显要好于 LSQ 的结果. 对序列 {x s t j ), y s t j ), j 1, 2,, k} 在时间域上作平均, 可以获得如下关系式, x 2 s + y 2 s a2 s 2 x s 2 + ȳ s 2 9e2 s 4 a2 s 1 + cos 2 i s ) + O e 2 ), 23) [ 1 sin 2 ω s sin 2 i s ] + O e 3 ). 24)

4 期许伟维等 : 天体测量法探测系外行星 433 假定已知行星质量 m p 恒星的平近点角 M 恒星的升交点经度 Ω 和恒星的近星点经度 ω, 可以由上述 3 个约束关系式 21) 23) 式, 近似给出恒星的轨道半长径偏心率和轨道倾角, 这样将对 7 维参数空间扫描问题降低到 4 维参数空间扫描, 可以显著提高计算效率. 当然, 在由观测时间序列计算上述平均值时, 会产生计算误差, 将会直接影响到上述约束条件的精度. 7 多行星系统的计算坐标速度法可以作为天体测量法探测系外行星的基本方法, 不仅是因为它可以用于小偏心率轨道计算, 而且较容易将其推广到多行星系统, 求解各个行星的参数位置速度及质量 ). 下面仅推导由 1 颗恒星和 2 颗行星组成的三体系统中恒星观测方程计算具体关系式, n 体问题的推导思路类似. 建立如图 4 所示惯性坐标系, 恒星与两颗行星质量分别为 M s m 1 和 m 2, O 为系统的质量中心, 空间位置矢量如图所示. 图 4 三体运动坐标示意图 Fig. 4 The coordinate s sketch of the three-body motion 3 个天体对应的运动方程为 GM s m i r i r s ) M s rs r i r s 3 i1 m i ri GM sm i r i r s ) r i r j 3 + Gm im j r i r j ) r i r j 3 i, j 1, 2 且 i j), 25) 因系统质心为坐标原点, 由质心运动定理知 : M s r s + m 1 r 1 + m 2 r 2 0, 26) 其中, r s x s, y s, z s ). 由行星引力导致恒星运动的观测方程可写为 M s δ 0 + ϱ δ ϖ + µ δ t t 0 ) ηt)) M s x s t) m 1 x 1 t) + m 2 x 2 t) M s α 0 + ϱ α ϖ + µ α t t 0 ) ξt)) M s y s t) m 1 y 1 t) + m 2 y 2 t). 27) 方程 25) 左边分别用 η t) ξ t) 表示, 可将其改写为 η t) m 1 x 1 t) + m 2 x 2 t) ξ t) m 1 y 1 t) + m 2 y 2 t), 28)

434 天文学报 57 卷上式中, x i t) y i t) 分别表示第 i 颗行星的位置坐标 x y 分量, 它们是运动方程初始条件 r 1 t 0 ) r 1 t 0 ) m 1 r 2 t 0 ) r 2 t 0 ) m 2 的函数. 仿照二体问题微分改正的做法, 假设初始估计值为 r 1 r 2 r 1 r 2 m 1 m 2, 上式在初始运动位置处展开, 对应观测方程的微分改正公式为 ηt) x i t) m i + i1 i1 j1 ξt) i1 j1 ) x i m i m m j j yi t) m i + i1 i1 j1 i, j 1, 2) m i y i m j ) i1 j1 m j 由 24) 25) 式可求出行星运动方程为 m i x i r j ) m i y i r j ) r j + r j + i1 j1 i1 j1 Ms + m i GM s r i + m ) j r j M s M s r i M s + m i r i + m 3 Gm j r i r j ) j r j r i r j 3 M s M s m i x i r j m i y i r j ) ) r j + r j +. 29) F i. 30) 令 ri v i, 则线性方程组 28) 的系数矩阵可由行星运动方程及其对应的状态转移矩阵联立求得, r i v i v i F i r i, r j, m i ) i j) ) ) d ri vi dt r j r ) j ) d ri vi dt v j v j ) d vi ) F i ri dt r j r i r j ) d vi ) F i ri dt v j r i v j ) ) d ri vi dt m j m j ) d vi ) F i + dt m j m j ) + ) + F i r i F ) ) i rj r j r j F ) ) i rj r j m j v j ) ) ri F + ) ) i rj r j m j, 31)

4 期许伟维等 : 天体测量法探测系外行星 435 状态转移方程涉及到的初始条件为 : 当 t t 0 时, r i H i j) r i, 0, r j 0 i j) v j v i v i H i j) 0,, r j v j 0 i j) r i m j v i m j 0 0 i, j 1, 2). 32) 这里, H 为 3 阶单位矩阵. 8 结论与讨论本文基于天体测量法给出的恒星空间位置参数序列, 讨论了在求解行星质量和轨道参数时涉及的动力学条件方程计算问题 ; 对于采用 Kepler 轨道根数和采用坐标速度描述形式分别给出了具体微分改正公式, 同时也进行了初步必要的仿真模拟计算 ; 通过天体力学二体问题运动关系式, 对恒星空间位置时间序列做平均分析, 推导出了待估参数间存在的 3 个约束条件关系式, 可以显著提升计算效率. 本文没有对恒星自行和视差等天体测量参数做解算研究, 以后若有实测数据, 可根据实际情况作具体处理 : 或单独求解, 或与行星轨道和质量一起求解. 参考文献 [1] Mayor M, Queloz D. Nature, 1995, 368: 355 [2] Michael P. The Exoplanet Handbook. New York: Cambridge University Press, 2011: 3-5 [3] Luisa G J, Luis A, Barbara P. MNRAS, 2014, 443: 260 [4] 姜明达, 肖东, 朱永田. 天文学进展, 2012, 30: 214 [5] Dou J P, Zhu Y T, Ren D Q. Chinese Journal of Nature, 2014, 36: 124 [6] 张牛, 季江徽. 天文学进展, 2009, 27: 14 [7] Lin D N C, Ida S. ApJ, 1997, 477: 781 [8] Papaloizou J C B, Terquem C. MNRAS, 2001, 325: 221 [9] 金文敬. 天文学进展, 2013, 31: 375 [10] Wright J T, Howard A W. ApJS, 2009, 182: 205 [11] Sara S. Exoplanets. Arizona: University of Arizona Press, 2010: 327-329 [12] Michael P, Joel H, Gaspar B, et al. ApJ, 2014, 797: 14 [13] Thiele T N. Astronomische Nachrichten, 1883, 104: 245

436 天文学报 57 卷 Exoplanet Detection by Astrometric Method XU Wei-wei 1,2,3,4 LIAO Xin-hao 1,3 ZHOU Yong-hong 1,3 XU Xue-qing 1,3 1 Shanghai Astronomical Observatory, Chinese Academy of Sciences, Shanghai 200030) 2 School of Physical Science and Technology, ShanghaiTech University, Shanghai 201210) 3 Key Laboratory of Planetary Sciences, Chinese Academy of Sciences, Shanghai 200030) 4 University of Chinese Academy of Sciences, Beijing 100049) ABSTRACT As we known, the exoplanets are mostly detected by the methods of radial velocity and transit, only one is found by the astrometric method. As the data of the gaia to be released, astrometry will become one of the most important method for detecting exoplanets gradually. Based on the position sequence of stars, this paper discusses the calculation of the equations of dynamics conditions involved in solving the mass and the orbit parameters of the planet. Due to the deficiency of the available theory orbital element method), we put forward a new method coordinate velocity method). The differential correction formulae of the two methods are presented, as well as the necessary simulation. In addition, the method established in this paper can be applied to the multi-planet system easily. Key words exoplanet, astrometry, celestial mechanics: orbit calculation and determination, methods: data analysis

4 期 许 伟 维 等 : 天 体 测 量 法 探 测 系 外 行 星 423 度 法 为 主 要 的 探 测 方 法 ; 2009 年 美 国 国 家 航 空 航 天 局 National Aeronautics and Space Administration, NASA) 建 造 的 开 普

4 期许伟维等 : 天体测量法探测系外行星 423 度法为主要的探测方法 ; 2009 年美国国家航空航天局 National Aeronautics and Space Administration, NASA) 建造的开普