CN109798921B - 一种星敏感器内方位元素室内校准方法 - Google Patents

Abstract

本发明公开了一种星敏感器内方位元素室内校准方法。该方法基于一个由气垫平台、一维单轴转台、单星星光模拟器、星敏感器、可旋转的星敏感器工装支架以及数据处理计算机组成的校准系统，其步骤包括：校准系统建模、数据采集以及数据处理。本发明采用了不依赖于外参数的算法模型，避免了因外部参数求解误差引入到内部参数的估计过程中，从而影响参数校准精度。本发明采用了一维单轴转台，降低了对校准设备的要求，同时规避了二维轴向转台使用过程中诸如星敏感器在内框中的安装误差、二维轴向转台两轴系不垂直引入的系统误差等问题。该方法减少了校准过程中复杂的星敏感器安装对准过程，使得整个校准过程变得相对简单。

Description

一种星敏感器内方位元素室内校准方法

技术领域

本发明属于航天测量技术领域，具体涉及一种星敏感器内方位元素室内校准方法。

背景技术

星敏感器是一种高精度的姿态敏感测量仪器，它以恒星为参照系、以星空为工作对象，利用其携带的星光相机拍摄星图，敏感多颗恒星给出多个参考矢量，通过星图预处理、星图识别、姿态解算等步骤解算出载体的姿态。星敏感器的校准是实现其姿态准确测量必不可缺少的重要环节。在星敏感器使用前，星敏感器的主点、焦距以及畸变等内部参数必须进行精确的测量，称为星敏感器校准。星敏感器校准一般分为地面校准和在轨校准两种方式，其中地面校准又包括外场校准和实验室校准，外场校准在大气通透、光污染小的地点对天顶晴朗夜空拍摄星图，提取星点坐标数据用于校准；实验室校准在星光实验室利用星光模拟器配合高精度二维轴向转台进行校准。其中实验室校准因为过程受控，是最基本也是精度最高的校准方式。实验室校准系统一般包括气垫平台、二维轴向转台、单星星光模拟器、星敏感器以及数据处理计算机。目前实验室校准系统存在的问题：

(a)二维轴向转台相对一维单轴转台制造成本更高，而且校准过程复杂；

(b)二维轴向转台使用过程中会引入诸如星敏感器在内框中的安装误差、二维轴向转台两轴系不垂直引入的系统误差等问题；

(c)即便在校准系统建模过程中将安装误差、测量系统误差等作为模型参数解算，但是此类校准方法使得外参数的求解误差反过来也影响星敏感器内部参数精度。

发明内容

本发明的目的是：针对上述校准方法存在的问题，提供一种基于单轴转台的星敏感器内方位元素室内校准方法，简化校准过程。在内部参数建模过程中，剥离外部参数的影响，采用非线性最小二乘法，通过线性化方程组迭代求解星敏感器光学系统内部参数。

本发明采用的技术方案是：一种星敏感器内方位元素室内校准方法，基于一个由气垫平台、一维单轴转台、单星星光模拟器、星敏感器、星敏感器工装支架(可旋转)以及数据处理计算机组成的校准系统(见图1)。星敏感器安装在工装支架内，工装支架固定在一维单轴转台上，而一维单轴转台和单星星光模拟器放置在气垫平台上。星敏感器视轴垂直于一维单轴转台的旋转轴，并且对准单星星光模拟器，星敏感器光学系统入瞳过一维单轴转台的旋转轴。以上安装无精度要求。校准的步骤如下：

1、校准系统建模

1.1坐标系

(a)一维单轴转台坐标系O_s-X_sY_sZ_s

将一维单轴转台的旋转轴定义为O_sZ_s轴，坐标系原点O_s在转台台面上，转台零位方向定义为O_sX_s轴，通过正交关系和右手系法则确定O_sY_s轴，参见图4。

(b)相机坐标系O’-X_oY_oZ_o

过星敏感器的光学系统中心O，并垂直于靶面Σ直线作为相机坐标系的O’Z_o轴，其中O’为垂线O’Z_o与靶面Σ的交点，O’X_o轴平行于靶面Σ的列方向，O’Y_o轴平行于靶面Σ的行方向，参见图2。

(c)靶面Σ坐标系O_f-X_fY_f

将探测器左上角坐标(0,0)作为靶面Σ坐标系原点O_f，O_fX_f轴平行于靶面Σ的列方向，O_fY_f轴平行于靶面Σ的行方向。当相机按逆时针方向旋转一角度时，记新的靶面Σ坐标系相对旋转前靶面Σ坐标系O_f-X_fY_f为O_f’-X_f’Y_f’，参见图3。

1.2针孔成像模型

在相机坐标系O’-X_oY_oZ_o下，星点S_i和S_j星光经过光学系统中心后，投影成像在星敏感器靶面星像点为S′_i和S′_j，对应在靶面Σ坐标系O_f-X_fY_f的坐标分别为P_i(x_i,y_i)和P_j(x_j,y_j)；假设星敏感器的焦距为f，主点O’在靶面Σ坐标系O_f-X_fY_f的坐标为(x_o,y_o)。光线S_iS_i′和S_jS_j′在一维单轴转台坐标系O_s-X_sY_sZ_s的单位方向为v_i和v_j。根据共线法则：

其中，R_so为相机坐标系O’-X_oY_oZ_o到一维单轴转台坐标系O_s-X_sY_sZ_s的旋转变换矩阵，

1.3畸变补偿模型

出于两个方面的考虑，一是，畸变补偿模型应该具备很强的泛化能力，不因不同星敏感器的畸变差异，导致模型不适用；二是，畸变分布应是一光滑曲面。采用常数c_x和c_y补偿星像点坐标P(x,y)的低阶畸变

通过焦距多项式修正高阶畸变，假设f(x,y)为星像点P(x,y)处的焦距，则：

式中P为多项式的阶数，a_st多项式系数。在畸变补偿后，式(2)被改写为：

其中焦距为：

2、数据采集

(a)将星敏感器安装在工装支架内，工装支架能够使得星敏感器绕光轴进行旋转，但对工装支架无精度要求；

(b)记录当前工装旋转位置，标记为第k组实验，假设共进行K＝180°/Δθ组实验，则1≤k≤K，其中0<Δθ≤45°为每次工装旋转角度；

(c)第k组实验数据采集：转动一维单轴转台，使单星星光模拟器的星像点扫过星敏感器靶面，每间隔Δα度为一个数据采集位置，在一个数据采集位置重复采集n次数据，n可取10～100。通过取平均得到星像点坐标P_ki(x_ki,y_ki)，并记录当时的转台转动读数α_ki，假设一共有N_k个采集位置，则1≤i≤N_k；

(d)将工装逆时针旋转一角度Δθ，进入下一组实验数据采集。

3、数据处理

假设R_so,k为第k组实验相机坐标系O’-X_oY_oZ_o到一维单轴转台坐标系O_s-X_sY_sZ_s的旋转变换矩阵，那么第k组实验数据有：

又因为v_ki＝(cosα_ki,sinα_ki,0)^T,v_kj＝(cosα_kj,sinα_kj,0)^T,则

于是得到一系列方程：

其中：

注意到式(8)中，常数(c_x,c_y)对主点位置已经进行了补偿，所以主点实际位置(x_o,y_o)在本发明中并不必要，仅作为一参考点。在本发明中,对主点位置无任何精度要求，不妨假设主点位置已经通过准直法获得。这样需要考虑校准的参数包括常值补偿量(c_x,c_y)和焦距多项式系数，共有

个参数，用参量X＝(c_x,c_y,a₀₀,a₁₀,…,a_P-1,P,a_PP)^T来表示。

下面采用迭代法求解方程组(6)：

(a)迭代过程初始化

预估计焦距

其中P_ki(x_ki,y_ki)和P_kj(x_kj,y_kj)为第k组数据靠近主点(x_o,y_o)的星像点，α_ki和α_kj为对应的转台角度读数。预估计畸变参数

和

(b)估计低阶畸变参数c_x、c_y和焦距参数a₀₀

设置多项式阶数P＝0，此时参量X＝(c_x,c_y,a₀₀)^T，将方程(7)在预估的参数附近线性展开：

写成矩阵形式：AΔX＝b,

其中

ΔX＝(Δc_x,Δc_y,Δa₀₀)^T，

令

采用最小二乘法，得到线性表达式(9)的解ΔX＝(A^TA)^-1A^Tb。为得到精确数值解，迭代求解焦距和畸变参数：

经过迭代5次后停止，固定畸变参数c_x和c_y，转入步骤(c)。

(c)检测迭代是否停止

将星间角间距误差作为算法评价指标：

如果Δr<T_Ang或者P≥T_P迭代停止，其中T_Ang为角间距误差阈值，与星敏感器设计参数相关，T_P为多项式阶数阈值，决定了畸变复杂度的上限；否则，执行步骤(d)。

(d)高阶畸变参数计算

焦距多项式阶数自加

为了防止算法出现过拟合，固定先前计算得到的主点位置和多项式系数，则求解的参数为X＝(a_0,P,a_1,P-1,…,a_P-1,P,a_PP),采用式(10)类似的方法进行求解，迭代5次以后，转入步骤(c)。

本发明与现有技术相比的优点是：

(a)采用了不依赖于外参数的系统模型，避免了因外部参数求解误差引入到内部参数的估计过程中，从而影响参数校准精度；

(b)采用了一维单轴转台而非二维轴向转台，降低了对标较设备的要求。规避了二维轴向转台使用过程中诸如星敏感器在内框中的安装误差、二维轴向转台两轴系不垂直引入的系统误差等问题；

(c)该方法减少了校准过程中复杂的星敏感器安装对准过程，使得整个校准过程变得相对简单。

附图说明

图1为本发明的一种星敏感器内方位元素室内校准方法涉及的校准系统；

图2为本发明的一种星敏感器内方位元素室内校准方法的相机坐标系；

图3为本发明的一种星敏感器内方位元素室内校准方法的靶面Σ坐标系；

图4为本发明的一种星敏感器内方位元素室内校准方法的一维单轴转台坐标系；

图5为本发明的一种星敏感器内方位元素室内校准方法的仿真的星点在像面的分布。

具体实施方式

下面对本发明的更多细节进行说明。

本发明提出一种基于一维单轴转台对星敏感器校准内部参数方法，利用旋转矩阵保持光线矢量张角不变的性质对测量方程进行建模。采用非线性最小二乘法，通过线性化方程组迭代求解星敏感器光学系统内部参数。这种方法采用了不依赖于外参数的系统模型，避免了因外部参数求解误差引入到内部参数的估计过程中，从而影响参数校准精度。同时采用了一维单轴转台，降低了对标较设备的要求，规避了二维轴向转台使用过程中诸如星敏感器在内框中的安装误差、二维轴向转台两轴系不垂直引入的系统误差等问题。该方法减少了校准过程中复杂的星敏感器安装对准过程，使得整个校准过程变得相对简单。下面详细说明本发明的步骤。

1、校准系统建模

1.1坐标系

(a)一维单轴转台坐标系O_s-X_sY_sZ_s

将一维单轴转台的旋转轴定义为O_sZ_s轴，坐标系原点O_s中心在转台台面上，转台零位方向定义为O_sX_s轴，通过正交关系和右手系法则确定O_sY_s轴，参见图4。

(b)相机坐标系O’-X_oY_oZ_o

过星敏感器的光学系统中心O，并垂直于靶面Σ直线作为相机坐标系的O’Z_o轴，其中O’为垂线O’Z_o与靶面Σ的交点，O’X_o轴平行于靶面Σ的列方向，O’Y_o轴平行于靶面Σ的行方向，参见图2。

(c)靶面∑坐标系O_f-X_fY_f

将探测器左上角坐标(0,0)作为靶面Σ坐标系原点O_f，O_fX_f轴平行于靶面Σ的列方向，O_fY_f轴平行于靶面Σ的行方向。当相机按逆时针方向旋转一角度时，记新的靶面Σ坐标系相对旋转前靶面Σ坐标系O_f-X_fY_f为O_f’-X_f’Y_f’，参见图3。

1.2针孔成像模型

在相机坐标系O’-X_oY_oZ_o下，星点S_i和S_j星光经过光学系统中心后，投影成像在星敏感器靶面星像点为S′_i和S′_j，对应在靶面Σ坐标系O_f-X_fY_f的坐标分别为P_i(x_i,y_i)和P_j(x_j,y_j)；假设星敏感器的焦距为f，主点O’在靶面Σ坐标系O_f-X_fY_f的坐标为(x_o,y_o)。光线S_iS′_i和S_jS′_j在一维单轴转台坐标系O_s-X_sY_sZ_s的单位方向为v_i和v_j。根据共线法则：

其中，R_so为相机坐标系O’-X_oY_oZ_o到一维单轴转台坐标系O_s-X_sY_sZ_s的旋转变换矩阵，

1.3畸变补偿模型

采用常数c_x和c_y补偿星像点坐标P(x,y)的低阶畸变

通过焦距多项式修正高阶畸变，假设f(x,y)为星像点P(x,y)处的焦距，则：

式中P为多项式的阶数，a_st多项式系数。在畸变补偿后，式(2)被改写为：

其中焦距为：

2、数据采集

(a)将星敏感器安装在工装支架内，工装支架能够使得星敏感器绕光轴进行旋转，但对工装支架无精度要求；

(b)记录当前工装旋转位置，标记为第k组实验，假设共进行K＝180°/45°＝4组实验，则1≤k≤4；

(c)第k组实验数据采集：转动一维单轴转台，使单星星光模拟器的星像点扫过星敏感器靶面，每间隔Δα＝1°为一个数据采集位置，在一个数据采集位置重复采集n＝10次数据。通过取平均得到星像点坐标P_ki(x_ki,y_ki)，并记录当时的转台转动读数α_ki，假设一共有N_k个采集位置，则1≤i≤N_k；

(d)将工装逆时针旋转一角度Δθ＝45°，进入下一组实验数据采集。

3、数据处理

假设R_so,k为第k组实验相机坐标系O’-X_oY_oZ_o到一维单轴转台坐标系O_s-X_sY_sZ_s的旋转变换矩阵，那么第k组实验数据有：

又因为v_ki＝(cosα_ki,sinα_ki,0)^T,v_kj＝(cosα_kj,sinα_kj,0)^T,则

于是得到一系列方程：

其中：

注意到式(8)中，常数(c_x,c_y)对主点位置已经进行了补偿，所以主点实际位置(x_o,y_o)在本发明中并不必要，仅作为一参考点。在本发明中,对主点位置无任何精度要求，不妨假设主点位置已经通过准直法获得。这样需要考虑校准的参数包括常值补偿量(c_x,c_y)和焦距多项式系数，共有

个参数，用参量X＝(c_x,c_y,a₀₀,a₁₀,…,a_P-1,P,a_PP)^T来表示。

下面采用迭代法求解方程组(6)：

(a)迭代过程初始化

预估计焦距

其中P_ki(x_ki,y_ki)和P_kj(x_kj,y_kj)为第k组数据靠近主点(x_o,y_o)的星像点，α_ki和α_kj为对应的转台角度读数。预估计畸变参数

和

(b)估计低阶畸变参数c_x、c_y和焦距参数a₀₀

设置多项式阶数P＝0，此时参量X＝(c_x,c_y,a₀₀)^T，将方程(7)在预估的参数附近线性展开：

写成矩阵形式：AΔX＝b,

其中

ΔX＝(Δc_x,Δc_y,Δa₀₀)^T，

令

采用最小二乘法，得到线性表达式(9)的解ΔX＝(A^TA)^-1A^Tb。为得到精确数值解，迭代求解焦距和畸变参数：

经过迭代5次后停止，固定畸变参数c_x和c_y，转入步骤(c)。

(c)检测迭代是否停止

将星间角间距误差作为算法评价指标：

如果Δr<T_Ang或者P≥T_P迭代停止，其中T_Ang为角间距误差阈值，与星敏感器设计参数相关，T_P为多项式阶数阈值，决定了畸变复杂度的上限；否则，执行步骤(d)。

(d)高阶畸变参数计算

焦距多项式阶数自加

为了防止算法出现过拟合，固定先前计算得到的主点位置和多项式系数，则求解的参数为X＝(a_0,P,a_1,P-1,…,a_P-1,P,a_PP),采用式(10)类似的方法进行求解，迭代5次以后，转入步骤(c)。

实施例：

本节通过仿真和实测数据验证方法有效性。具体的实施过程的范例见以上“具体实施方式”这一节。为验证方法的泛化能力，拟合使用数据(用于计算内部参数)和验证使用数据(用于计算星间角间距误差)不重叠。

在第一个测试用例中，仿真的星敏感器的基本参数为：

视场：Φ20°

像素阵列：1536×1536

主点位置：(758，775)

像元尺寸：5.5um×5.5um

焦距：24.0955mm

假设星点质心噪声为均值为零，标准差为0.05像素的高斯分布噪声，由于每个数据采集位置采集次数为n＝10，经过平均后噪声水平下降，即0.016像素。同时，假设Δx和Δy为星像点P(x,y)在靶面Σ坐标系的O_fX_f方向和O_fY_f方向的畸变误差，按畸变理论公式：

式中：

径向畸变参数：q₁＝7.5E-9，q₂＝-5.0E-15，q₃＝-2.0E-21；切向畸变参数：p₁＝7.5E-7，p₂＝-5.5E-7。通过仿真，得到星点在像面的分布参见图5。为了验证本发明方法的有效性，在实施例中采用星间角间距误差作为算法评价指标，见式(11)。焦距多项式阶数阈值T_P＝5，数据采集过程中设置Δθ＝45°，Δα＝1°。最终，得到仿真数据处理结果见表1和表2。

表1仿真数据处理结果的星间角间距误差

表2仿真数据计算得到的畸变系数

c<sub>x</sub>	c<sub>y</sub>	a<sub>00</sub>	a<sub>01</sub>	a<sub>10</sub>	a<sub>02</sub>
						-42.6	31.1	4380.6	1.10E-03	-2.06E-03	-3.31E-05
a<sub>11</sub>	a<sub>20</sub>	a<sub>03</sub>	a<sub>12</sub>	a<sub>21</sub>	a<sub>30</sub>
						4.51E-06	-3.63E-05	-7.79E-10	2.76E-09	7.27E-10	4.65E-09
a<sub>04</sub>	a<sub>13</sub>	a<sub>22</sub>	a<sub>31</sub>	a<sub>40</sub>	a<sub>05</sub>
						2.17E-11	-1.94E-10	5.30E-11	1.65E-10	3.93E-11	-7.90E-17
a<sub>14</sub>	a<sub>23</sub>	a<sub>32</sub>	a<sub>41</sub>	a<sub>50</sub>
						-1.30E-11	2.76E-11	1.30E-11	-2.76E-11	-5.39E-15

在第二个测试用例中，将一星敏放在一维单轴转台上进行测试，设置标较参数Δθ＝45°，Δα＝0.2°,T_P＝5。星敏感器的基本参数为：

视场：Φ20°

像素阵列：1536×1536

主点位置：(765，745)(准直法实测)

像元尺寸：5.5um×5.5um

焦距：24.0000mm(光学系统设计理论值)

单轴转台定位精度：0.2”，大约相当于0.005个像素，像元定位误差(经过多次采集平均)的定位误差约0.05像素。

表3实测数据处理结果的星间角间距误差

表4实测数据计算得到的畸变系数

c<sub>x</sub>	c<sub>y</sub>	a<sub>00</sub>	a<sub>01</sub>	a<sub>10</sub>	a<sub>02</sub>
						-28.4	-12.8	4377.0	-3.69E-04	6.43E-04	-2.09E-06
a<sub>11</sub>	a<sub>20</sub>	a<sub>03</sub>	a<sub>12</sub>	a<sub>21</sub>	a<sub>30</sub>
						-3.20E-06	-2.77E-06	5.94E-09	5.07E-08	7.63E-09	-1.53E-08
a<sub>04</sub>	a<sub>13</sub>	a<sub>22</sub>	a<sub>31</sub>	a<sub>40</sub>	a<sub>05</sub>
						-7.36E-12	-8.76E-10	-3.65E-10	7.34E-10	7.55E-11	7.84E-15
a<sub>14</sub>	a<sub>23</sub>	a<sub>32</sub>	a<sub>41</sub>	a<sub>50</sub>
						2.80E-12	3.69E-12	-2.09E-12	-3.07E-12	1.07E-14

仿真数据处理结果显示星间角间距误差(1σ)为0.65”(采用4阶焦距多项式)，相当于0.014像素，与仿真质心算法噪声水平相一致。实测数据处理结果显示星间角间距误差(1σ)为2”(采用4阶焦距多项式)，相当于0.043像素，与实测质心算法噪声水平0.05个像元相一致。仿真数据和实测数据均表明本发明校准精度完全满足使用要求。

Claims (1)

1.一种星敏感器内方位元素室内校准方法，该方法基于一个由气垫平台、一维单轴转台、单星星光模拟器、星敏感器、可旋转的星敏感器工装支架以及数据处理计算机组成的校准系统，星敏感器安装在工装支架内，工装支架固定在一维单轴转台上，而一维单轴转台和单星星光模拟器放置在气垫平台上，星敏感器视轴垂直于一维单轴转台的旋转轴，并且对准单星星光模拟器，星敏感器光学系统入瞳过一维单轴转台的旋转轴，其特征在于：校准步骤如下：

1.1校准系统建模

1.1.1坐标系

(a)一维单轴转台坐标系O_s-X_sY_sZ_s

将一维单轴转台的旋转轴定义为O_sZ_s轴，坐标系原点O_s在转台台面上，转台零位方向定义为O_sX_s轴，通过正交关系和右手系法则确定O_sY_s轴；

(b)相机坐标系O’-X_oY_oZ_o

过星敏感器的光学系统中心O，并垂直于靶面∑的直线作为相机坐标系O’Z_o轴，其中O’为垂线O’Z_o与靶面∑的交点，O’X_o轴平行于靶面∑的列方向，O’Y_o轴平行于靶面∑的行方向；

(c)靶面∑坐标系O_f-X_fY_f

将探测器左上角坐标(0，0)作为靶面∑坐标系原点O_f，O_fX_f轴平行于靶面∑的列方向，O_fY_f轴平行于靶面∑的行方向；

1.1.2针孔成像模型

在相机坐标系O’-X_oY_oZ_o下，星点S_i和S_j星光经过光学系统中心后，投影成像在星敏感器靶面星像点为S′_iS′_j，对应在靶面∑坐标系O_f-X_fY_f的坐标分别为P_i(x_i，y_i)和P_j(x_j，y_j)，假设星敏感器的焦距为f，主点O’在靶面∑坐标系O_f-X_fY_f的坐标为(x_o，y_o)，光线S_iS′_i和S_jS′_j在一维单轴转台坐标系O_s-X_sY_sZ_s的单位方向为v_i和v_j，根据共线法则：

其中，R_so为相机坐标系O’-X_oY_oZ_o到一维单轴转台坐标系O_s-X_sY_sZ_s的旋转变换矩阵，

1.1.3畸变补偿模型

采用常数c_x和c_y补偿星像点坐标P(x，y)的低阶畸变

通过焦距多项式修正高阶畸变，假设f(x，y)为星像点P(x，y)处的焦距，则：

式中P为多项式的阶数，a_st多项式系数，在畸变补偿后，式(2)被改写为：

其中焦距为：

1.2数据采集

(a)将星敏感器安装在工装支架内，工装支架能够使得星敏感器绕光轴进行旋转，但对工装支架无精度要求；

(b)记录当前工装旋转位置，标记为第k组实验，假设共进行K＝180°/Δθ组实验，则1≤k≤K，其中0＜Δθ≤45°为每次工装旋转角度；

(c)第k组实验数据采集：转动一维单轴转台，使单星星光模拟器的星像点扫过星敏感器靶面，每间隔Δα度为一个数据采集位置，在一个数据采集位置重复采集n次数据，n可取10～100，通过取平均得到星像点坐标P_ki(x_ki，y_ki)，并记录当时的转台转动读数α_ki，假设一共有N_k个采集位置，则1≤i≤N_k；

(d)将工装逆时针旋转一角度Δθ，进入下一组实验数据采集；

1.3数据处理

假设R_so，k为第k组实验相机坐标系O’-X_oY_oZ_o到一维单轴转台坐标系O_s-X_sY_sZ_s的旋转变换矩阵，那么第k组实验数据有：

又因为v_ki＝(cosα_ki，sinα_ki，0)^T，v_kj＝(cosα_kj，sinα_kj，0)^T，则

于是得到一系列方程：

其中：

注意到式(8)中，常数(c_x，c_y)对主点位置已经进行了补偿，所以主点实际位置(x_o，y_o)并不必要，仅作为一参考点，对主点位置无任何精度要求，不妨假设主点位置已经通过准直法获得，这样需要考虑校准的参数包括常值补偿量(c_x，c_y)和焦距多项式系数，共有

个参数，用参量X＝(c_x，c_y，a₀₀，a₁₀，...，a_P-1，P，a_PP)^T来表示；

下面采用迭代法求解方程组(6)：

(a)迭代过程初始化

预估计焦距

其中P_ki(x_ki，y_ki)和P_kj(x_kj，y_kj)为第k组数据靠近主点(x_o，y_o)的星像点，α_ki和α_kj为对应的转台角度读数，预估计畸变参数

和

(b)估计低阶畸变参数c_x、c_y和焦距参数a₀₀

设置多项式阶数P＝0，此时参量X＝(c_x，c_y，a₀₀)^T，将方程(7)在预估的参数附近线性展开：

写成矩阵形式：AΔX＝b，

其中

ΔX＝(Δc_x，Δc_y，Δa₀₀)^T，

；

令

采用最小二乘法，得到线性表达式(9)的解ΔX＝(A^TA)^-1A^Tb，为得到精确数值解，迭代求解焦距和畸变参数：

经过迭代5次后停止，固定畸变参数c_x和c_y，转入步骤(c)；

(c)检测迭代是否停止

将星间角间距误差作为算法评价指标：

如果Δr＜T_Ang或者P≥T_P迭代停止，其中T_Ang为角间距误差阈值，与星敏感器设计参数相关，T_P为多项式阶数阈值，决定了畸变复杂度的上限；否则，执行步骤(d)；

(d)高阶畸变参数计算

焦距多项式阶数自加

为了防止算法出现过拟合，固定先前计算得到的主点位置和多项式系数，则求解的参数为X＝(a_0，P，a_1，P-1，...，a_P-1，P，a_PP)，采用式(10)类似的方法进行求解，迭代5次以后，转入步骤(c)。

Images