面向天体测量的矮行星Haumea相对测光试验中，告诉了我们些什么？

　　太阳系柯伊伯带天体轨道通常在海王星轨道以外，由于其日心距离远，温度低，它们是行星吸积盘原始遗迹之一，因此携带了有关于行星形成的物理和化学的宝贵信息。目前已知的柯伊伯带天体超过4000颗。
　　2004-2005年，Brown团队和Ortiz团队声称发现了矮行星Haumea。它被临时命名为2003 EL61，并于2006年9月7日由官方小行星收录，编号为136108。
　　Haumea是最大的柯伊伯带天体之一，并且是继冥王星和 (136472) Makemake后最亮的海外天体，其轨道周期为284地球年，近日点35AU，远日点51AU，绝对星等0.358 mag，半长轴43.3 AU，偏心率0.19，轨道倾角28 28 。
　　光变曲线是获取天体物理性质的主要来源，为小天体和行星的起源和演化的研究提供了重要线索。例如，导出Haumea自转光变曲线的周期为3.9154   0.0002 hr，峰峰值幅度为0.28   0.02 mag，颜色B-V = 0.97   0.03 mag。
　　使用夏威夷大学2.2m望远镜观测，观测时使用了Kron-Cousins测光系统的B，V，R和I滤光片，导出周期为3.9155   0.0001 hr的双峰自转光变曲线，其峰峰值幅度为0.29   0.02 mag。因自转周期短，Haumea呈旋转椭球体 (Jacobi ellipsoid)，假设其已达到流体静力平衡，结合这两个特性推测其质量密度ρρ~2500 kg m-3。
　　观测资料和数据处理
　　我们使用紫金山天文台姚安观测站80 cm望远镜（详细参数如表1）对矮行星Haumea进行两轮测光观测，观测时使用Johnson-Cousins I滤光片，具体的观测资料情况如表2。
　　利用MaxIm DL软件对获得的原始图像做预处理，主要包括本底、暗场和平场校正。由于反射光在CCD硅片薄层间来回反射，形成干涉条纹，导致图像背景的不均匀。
　　尤其对于Johnson-Cousins I滤光片，CCD硅片的吸收率随着波长的增加迅速下降，有相当一部分的低能光子在没有完全被吸收之前，在CCD硅片的薄层来回反射好几次，产生干涉条纹，不利于精确的光度测量。因此需要对拍摄的图像进行处理，剔除干涉条纹，具体的步骤如下。
　　首先，从每晚观测的CCD图像中选择1或2幅图像。接着，对筛选的CCD图像（剔除了3σ3σ）使用＂3-2＂公式（3倍中值- 2倍均值）进行天空背景的估计，得到每一幅图像各自的背景估计值。
　　所有筛选的CCD图像减去各自的背景估计值后，将这些图像每一个对应的灰度值取中值运算，将最后的中值结果作为干涉条纹模板对应像素的灰度值。
　　对于构造的干涉条纹模板，需要乘上一个比例因子以获得归一化干涉条纹模板，该比例因子为360 s除以待去除干涉条纹图像的曝光时间。最后，从待去除干涉条纹的原始图像中减去归一化干涉条纹模板，得到剔除了干涉条纹的图像。典型的处理结果如图1。
　　光度测量数据
　　本文使用相对测光技术获取目标的相对仪器星等。由于两轮观测时视场中的星象较少，而且在每轮观测期间，目标仅运动了一小段距离，因此我们可以使用相同的参考星进行测光。
　　具体地，在第1轮观测时，选择5颗参考星，在第2轮观测时，选择3颗参考星。所有的参考星都是恒星且光度测量稳定。
　　两轮观测中典型的目标和参考星的图像如图2，其中，目标Haumea使用绿色圆圈标记，参考星使用红色圆圈标记。
　　我们使用二维高斯定心算法获取目标星和参考星的像素坐标。接着，采用＂3-2＂公式对目标星和参考星周围环形区域天空背景值进行估计。
　　使用photutils测光库提供的孔径测光程序对目标星和参考星进行光度测量。photutils孔径测光程序允许我们使用一系列不同的孔径，选择光度测量精度最高的孔径作为最优孔径，得到每个夜晚目标星相对仪器星等。
　　误差主要来源于孔径的目标信噪比，同时我们也考虑了目标与参考星之间因颜色不同造成的误差。
　　由于没有进行标准测光，为了将两轮观测的测光结果联合起来进行周期分析，我们分别拟合两轮观测求出的目标相对仪器星等，拟合结果之间存在一个常数偏差，最后使用这个常数改正求出的目标相对仪器星等。每次拟合的是一轮多晚观测的数据，而不是每个夜晚的观测数据。
　　自转周期
　　为了测定准确的Haumea自转周期，本文使用两种周期求解方法，分别是相位弥散最小化方法以及Lomb-Scargle 周期图法。两种方法求出目标Haumea自转周期结果和讨论分析如下。
　　相位弥散最小化方法由Stellingwerf[9]在1978年提出。相位弥散最小化方法通过求统计量θθ的最小值来获取周期，θθ定义为：
　　θ=s2σ2θ=s2σ2(1)
　　σ2= Ni=1(xi x )2N 1，s2= Mj=1(nj 1)s2j Mj=1(nj M)σ2= i=1N(xi x )2N 1，s2= j=1M(nj 1)sj2 j=1M(nj M)(2)
　　(1)式中σ2σ2为观测的N个离散数据点的方差，(2)式中xixi是指在titi时刻观测的目标星等，x x 为N个数据点的均值。从N个数据点中选择M个有着相同相位ΦiΦi的样本，每个样本中包含数据点的个数为njnj，s2jsj2为第j个样本的方差。
　　通常，θθ值在0到1之间。如果求出的是正确的周期，那么θθ会达到一个局部最小值，接近0。使用相位弥散最小化方法计算得到的周期图如图3所示。在频率为12.25939 (cycles/day) 处达到局部最小值，最小值约为0.17。根据自转频率ΩΩ与自转周期PΡ之间的关系Ω=24/PΩ=24/Ρ，可以求出自转周期PΡ=1.95768 hr。
　　自转光变曲线分析
　　对于像Haumea这种因快速旋转而变形拉长成椭球体的TNO (Trans-Neptunian Object)，其完整的光变曲线应包含两个极大值和两个极小值，极大值代表天体的前后两个侧视面，而最小值是两端的视图。因此，将2.2节求解的Haumea自转周期拉长至两倍，求得Haumea的自转周期T=3.9154   0.0002 hr。
　　我们使用一般的线性最小二乘法拟合测光数据（剔除了3σ3σ），拟合模型为
　　mi=c+ j=13[ajsin(j2πTti)+bjcos(j2πTti)]mi=c+ j=13[ajsin(j2πTti)+bjcos(j2πTti)]， (5)
　　其中，titi为观测时刻的JD时间（校正了光行时）。各系数拟合结果和误差如表3。观测数据点和拟合模型光变曲线如图5，图中显示了两个极大值分别为0.12 mag和0.10 mag，两个极小值分别为0.14 mag和0.09 mag。对应的χ2χ2=252，拟合结果服从自由度为205的χ2χ2分布。χ2χ2定义为
　　χ2= i=1N[yi miσi]2χ2= i=1N[yi miσi]2， (6)
　　其中，N为拟合数据点总个数；yiyi为第ii个数据点的相对仪器星等；mimi为(5)式计算的结果；σiσi为第ii个数据点的测量误差。
　　光变对位置测量影响的讨论
　　对于像Haumea这种光度变化幅度较大的天体，光度会对位置测量造成一定的影响。由于Namaka较暗，现只考虑Haumea与卫星Hi’iaka构成的双星系统（如图6）。在天空平面上，设主星与卫星之间的角距离为dd；质心与光心之间的角距离为ΔΔ；P和S分别代表主星和卫星；CMCM和CLCL分别是质心和光心。
　　根据杠杆原理，
　　Δ=[11+rL 11+rM]dΔ=[11+rL 11+rM]d， (7)
　　其中，rL=Lp/LsrL=Lp/Ls；rM=mp/msrM=mp/ms；LpLp和LsLs分别是主星和卫星的亮度；mpmp和msms分别是主星和卫星的质量。假设卫星Hi’iaka的光变幅度不变，ΔΔ的变化主要由Haumea光变引起，可以根据误差估计公式估计ΔΔ受光变影响的变化范围。令
　　k=11+rL 11+rMk=11+rL 11+rM， (8)
　　则有Δ=kdΔ=kd，那么ΔΔ受光度的变化影响：
　　Δ′=(k+δk)dΔ＂=(k+δk)d， (9)
　　kk的变化只与主星和卫星的亮度变化有关，因此有
　　δk LsLΣδLΣLΣδk LsLΣδLΣLΣ. (10)
　　其中，LΣLΣ是主星与卫星的亮度总和，卫星Hi’iaka的亮度是主星Hamuea的5.9%，因此rL=16.9rL=16.9，Ls/LΣ=1/17.9Ls/LΣ=1/17.9。
　　rM=223.80rM=223.80。而根据星等转换关系式有δLΣ/LΣ=10 0.4 Δm 1δLΣ/LΣ=10 0.4 Δm 1，其中，ΔmΔm为Haumea的光变曲线上对应的星等值。
　　从IMCCE历表中可以获得任意时刻主星与卫星之间的距离dd。在整个观测期间，选取Haumea与卫星Hi’iaka的距离最远（d 1.446d 1.446′′）的时间节点，每隔1min取数据点，一共300个数据点，覆盖了Haumea的一个自转周期（T = 3.9154 hr）。
　　根据2.3节求得的Haumea光变曲线图（图5），我们可以获得Haumea光度的变化对位置测量的影响曲线图（如图7）。其中，横坐标为Haumea的自转相位，与图5的横坐标相对应，纵坐标为该时刻因光度变化，Haumea位置测量所受到的影响。可以看出，Haumea光度的变化对位置测量的影响最大可达到-9 ~ 9 mas。因此，利用Haumea的光变特征改正其位置测量是有必要的。
　　总结与展望
　　我们对紫金山天文台姚安观测站80cm望远镜拍摄的矮行星Haumea CCD图像剔除了干涉条纹，使用相对测光技术获得Haumea的相对仪器星等，通过联合两轮观测的资料探究了其光变情况。
　　使用相位弥散最小化方法以及Lomb-Scargle周期图法计算Haumea的自转周期，结果显示了一条明显的双峰自转光变曲线，其周期为3.9154   0.0002 hr，峰峰值幅度为0.26   0.01 mag。最终估计了Haumea由于自转产生的光度变化对位置测量的影响最大为-9 mas ~ 9 mas。在后续的工作中，我们将结合Haumea光度测量与位置测量，利用其光度变化特征改进位置测量精度。参考文献
　　[1] ORTIZ J L, SANTOS-SANZ P, SICARDY B, et al. The size, shape, density and ring of the dwarf planet Haumea from a stellar occultation[J]. Nature, 2017, 550:219-223.
　　[2] RABINOWITZ D L, BARKUME K, BROWN M E, et al. Photometric observations constraining the size, shape, and albedo of 2003 EL61, a rapidly rotating, Pluto-sized object in the Kuiper Belt[J]. The Astrophysical Journal, 2006, 639:1238-1251.
　　[3] LACERDA P, JEWITT D C, PEIXINHO N. High-precision photometry of extreme KBO 2003 EL61[J]. The Astronomical Journal, 2008, 135:1749-1756.
　　[4] BROWN M E, VAN DAM M A, BOUCHEZ A H, et al. Satellites of the largest Kuiper Belt objects[J]. The Astrophysical Journal, 2006, 639(1):L43-L46.
　　[5] TORRES G. Orbital solution and dynamical masses for the nearby binary system GJ67AB[J]. Monthly Notices of the Royal Astronomical Society, 2022, 514:5183-5191.
　　[6] 郭碧峰, 彭青玉, 林孚荣. Astrometrica与MaxIm DL在天体测量的应用与技术研究[J]. 天文研究与技术, 2021, 18(2):195-202.
　　[7] 郭碧峰, 彭青玉, 尚一佳, 等. CCD图像相对测光中的技术研究与精度回报现象[J]. 天文研究与技术, 2022, 19(2):111-117.
　　[8] THIROUIN A. Study of trans-neptunian objects using photometric techniques and numerical simulations[D]. Granada:Universidaa de Alicante, 2013.
　　[9] STELLINGWERF R F. Period determination using phase dispersion minimization[J]. The Astrophysical Journal, 1978, 224:953-960.
　　[10] LOMB N R. Leastsquares frequency analysis of unequally spaced data[J]. Astrophysics and Space Science, 1976, 39:447-462.
　　[11] SCARGLE J D. Studies in astronomical time series analysis. II. statistical aspects of spectral analysis of unevenly spaced data[J]. The Astrophysical Journal, 1982, 263:835-853.