CN100376883C - 一种基于像素频率的星敏感器校准方法 - Google Patents

Abstract

本发明航天测量技术，涉及对星敏感器校准方法的改进。本发明使用由气垫平台、单星星光模拟器、星敏感器、2维轴向转台以及数据处理计算机组成的校准系统，其步骤是：建立点扩散函数模型，建立质心算法公式，建立质心算法像素频率误差模型，数据采集。本发明提出了更为实用的点扩散模型；提高了星敏感器的校准精度；校准模型简单，利于计算。

Description

一种基于像素频率的星敏感器校准方法

技术领域

本发明航天测量技术，涉及对星敏感器校准方法的改进。

背景技术

星敏感器是一种利用恒星观测，为空间飞行器提供高精度姿态信息的航天测量仪器。在星敏感器投入使用前，其主点、焦距以及畸变系数等内部参数必须进行精确的测量，称为星敏感器校准。常见的星敏感器地面校准方法主要有两种：一种是在星光实验室利用星光模拟器配合高精度2轴转台进行数据采集和校准；另一种是在大气环境良好的地点对晴朗夜空进行拍摄以获得数据并进行校准。如图1所示为实验室内部星敏感器校准系统示意图，该校准系统主要有气垫平台，单星星光模拟器，星敏感器，2维轴向转台以及数据处理计算机组成。目前，实验室内部的星敏感器校准方法为精度最高的校准方法，通常的星敏感器校准精度可以达到1/10～1/20像素。这一校准水平认为是达到了星敏感器星点成像质心算法的误差水平，但是根据对星敏感器质心算法误差分析表明，该误差并不是简单的随机误差，而是由算法本身的系统误差和随机误差组成，其系统误差是一个以像素为变化频率的近似正弦偏差。目前的星敏感器校准过程中往往忽略了质心算法的系统误差问题，而把它简单的看成是随机噪声降低了工作精度。

发明内容

本发明的目的是：提出一种对星敏感器质心算法系统误差进行建模和校准的方法。该方法在星敏感器校准系统的内方元素校准完成后，能进一步提高星敏感器的校准精度。由于星敏感器质心算法系统误差是以像素频率出现的，因此称该方法为像素频率的高精度星敏感器校准。经过像素频率校准，可以将星敏感器星点的质心精度从1/10～1/20提高到1/100～1/200像素精度，达到星点质心精度极限水平。

本发明的技术方案是：一种基于像素频率的星敏感器校准方法，使用由气垫平台、单星星光模拟器、星敏感器、2维轴向转台以及数据处理计算机组成的校准系统，其特征在于，

1、星敏感器质心算法；

1.1、建立点扩散函数模型，星像在图像传感器靶面上的分布称为点扩散函数，通过计算弥散圆斑的质心来获得亚像素精度的星点位置坐标，点扩散函数模型见式(1)：

$I(x,y)=\left\{\begin{array}{cc}{I}_0\exp (-(x^2+y^2))& x^2+y^2\&GreaterEqual;r^2\\ I_0& x^2+y^2<r^2\end{array}\right.---\left(1\right)$

式中，I₀为中心信号输出，x，y为像素坐标，r为截断半径，该模型是点扩散函数在某个高度截断后的形状；

1.2、建立质心算法公式为：

$x_c=\frac{\&Integral;{\&Integral;}_A\mathrm{xI}(x,y)\mathrm{dxdy}}{\&Integral;{\&Integral;}_{A}I(x,y)\mathrm{dxdy}},$ $y_c=\frac{\&Integral;{\&Integral;}_A\mathrm{yI}(x,y)\mathrm{dxdy}}{\&Integral;{\&Integral;}_{A}I(x,y)\mathrm{dxdy}}---\left(2\right)$

式中，x_c，y_c为星点辐射中心，A表示大于阈值的区域，x，y为靶面像素坐标，I(x，y)为辐射分布函数，经过靶面像素的离散化后，质心计算公式为：

${\tilde{x}}_c=\frac{{\mathrm{\&Sigma;}}_{k=1}^n{x}_k{I}_k}{{\mathrm{\&Sigma;}}_{k=1}^n{I}_k},$ ${\tilde{y}}_c=\frac{{\mathrm{\&Sigma;}}_{k=1}^n{y}_k{I}_k}{{\mathrm{\&Sigma;}}_{k=1}^n{I}_k}---\left(3\right)$

式中，

为离散化后的星点中心，n为阈值区域内像素个数，k为像素序号，x_k，y_k第k个像素的坐标，I_k为第k个像素内积累的光电子对应亮度输出；

2、建立质心算法像素频率误差模型；

2.1、质心算法像素频率误差模型见式(4)：

A_x(sin(2πx_p+2πB_x)+A_y(sin(2πy_p+2πB_y)-sin(2πB_x)-sin(2πB_y))  (4)

式中，Ax、Ay为能量误差系数，Bx，By为相位误差系数，(x_p，y_p)为一个像素内的位置坐标；在整个靶面上，Ax，Ay为常数，Bx，By为2阶变化；

2.2、对上述系数在整个靶面上进行建模为：

A_x＝A_y＝T

$B_x=\left(\frac{x}{N}\right)b_1-{\left(\frac{x}{N}\right)}^2{b}_2---\left(5\right)$

$B_y={\left(\frac{y}{N}\right)}^2{b}_1-{\left(\frac{y}{N}\right)}^2{b}_2$

式中，T为常数，N为靶面行或者列像素数，(x，y)为靶面像素整数坐标，b1，b2为相位偏差系数；

2.3、根据公式(4)得到像素频率偏差函数：

f(T，Bx，By)＝

T(sin(2πx_p+2πB_x)+T(sin(2πy_p+2πB_y)-sin(2πB_x)-sin(2πB_y))  (6)

$\mathrm{\&Delta;d}=\left(\begin{array}{ccc}\frac{\&PartialD;f}{\&PartialD;T}& \frac{\&PartialD;f}{\&PartialD;\mathrm{Bx}}& \frac{\&PartialD;f}{\&PartialD;\mathrm{By}}\end{array}\right)\left(\begin{array}{c}\mathrm{\&Delta;T}\\ \mathrm{\&Delta;Bx}\\ \mathrm{\&Delta;By}\end{array}\right)---\left(7\right)$

式中，Δd为参数估计误差，ΔT为能量误差，ΔBx、ΔBy分别为x和y方向相位偏差，采用最小二乘迭代估计方法，可以求得以上3个参数的估计值；然后计算得到各个像素Bx和By的估计值，并计算整个靶面的相位偏差系数b1和b2；根据公式(5)有：

$\mathrm{\&Delta;Bx}=\left[\begin{array}{cc}\frac{x}{N},& -{\left(\frac{x}{N}\right)}^2\end{array}\right]\left[\begin{array}{c}\mathrm{\&Delta;}{b}_1\\ \mathrm{\&Delta;}{b}_2\end{array}\right]---\left(8\right)$

$\mathrm{\&Delta;By}=\left[\begin{array}{cc}\frac{y}{N},& -{\left(\frac{y}{N}\right)}^2\end{array}\right]\left[\begin{array}{c}\mathrm{\&Delta;}{b}_1\\ \mathrm{\&Delta;}{b}_2\end{array}\right]$

然后，联立两方程采用最小二乘法估计出b1和b2；

3、数据采集；

3.1、在一个像素内的采集，采集点阵为9×9点阵；根据9×9的点阵数据计算位置偏差，然后来提取该像素的模型参数；

3.2、选择整个靶面的像素进行9×9采样，所采用的星敏感器靶面为1024×1204像素阵列，x方向采样位置为：1、129、257、385、513、641、769、897、1024，采样间距为128；在y方向采样间距也为128，得到9×9的靶面采集像素位置；将这些位置上的测量值和估计值之间的偏差代入公式(8)中计算得到b1和b2。

本发明的优点是：

1、提出了更为实用的点扩散模型；

2、提高了星敏感器的校准精度；

3、校准模型简单，利于计算。

附图说明

图1是星敏感器校准系统组成示意图。

图2是星敏感器点扩散函数示意图。

图3是星敏感器像素频率误差示意图。

图4是像素内数据采集示意图。

具体实施方式

下面对本发明做进一步详细说明。

星敏感器质心算法。

星敏感器是通过对恒星进行拍摄和识别来完成姿态计算的，由于普通恒星的张角远远小于一个像素的对应角，因此聚焦成像的结果将是成像在在一个像素内。为了提高星敏感器质心定位精度到亚像素水平，通常采用离焦技术使得星点在图像传感器靶面上形成一个弥散圆斑，星像在图像传感器靶面上的分布称为点扩散函数。如图2所示为X方向点扩散示意图，I表示辐射能量分布。通过计算弥散圆斑的质心来获得亚像素精度的星点位置坐标。通过分析恒星的离焦和衍射作用，这里提出一种简单的点扩散函数模型：

$I(x,y)=\left\{\begin{array}{cc}{I}_0\exp (-(x^2+y^2))& x^2+y^2\&GreaterEqual;r^2\\ I_0& x^2+y^2<r^2\end{array}---\left(1\right)\right.$

这里，I0为中心信号输出，x，y为像素坐标，r为截断半径。该模型是点扩散函数在某个高度截断后的形状。

常用的质心算法公式为：

$x_c=\frac{\&Integral;{\&Integral;}_A\mathrm{xI}(x,y)\mathrm{dxdy}}{\&Integral;{\&Integral;}_{A}I(x,y)\mathrm{dxdy}},$ $y_c=\frac{\&Integral;{\&Integral;}_A\mathrm{yI}(x,y)\mathrm{dxdy}}{\&Integral;{\&Integral;}_{A}I(x,y)\mathrm{dxdy}}---\left(2\right)$

这里，x_c，y_c为星点辐射中心，A表示大于阈值的区域，x，y为靶面像素坐标，I(x，y)为辐射分布函数。经过靶面像素的离散化后，质心计算公式为：

${\tilde{x}}_c=\frac{{\mathrm{\&Sigma;}}_{k=1}^n{x}_k{I}_k}{{\mathrm{\&Sigma;}}_{k=1}^n{I}_k},$ ${\tilde{y}}_c=\frac{{\mathrm{\&Sigma;}}_{k=1}^n{y}_k{I}_k}{{\mathrm{\&Sigma;}}_{k=1}^n{I}_k}---\left(3\right)$

这里，

为离散化后的星点中心，n为阈值区域内像素个数，k为像素序号，x_k，y_k第k个像素的坐标，I_k为第k个像素内积累的光电子对应亮度输出。

影响质心算法精度的因素主要有两个，一个是镜头不完善引起的弥散圆的变形；另一个则是由于像素对辐射分布的离散化作用，用像素的几何中心代替了辐射分布的平衡点。如图3所示为靶面不同位置采集得到像素频率误差示意图。由于x和y方向对质心算法而言是独立的，因此这里仅采用y_p＝0，x_p方向的误差表示。图3中的4个小图，左上，右上，左下，右下分别对应于像素(205，512)，(410，512)，(615，512)和(820，512)。

质心算法像素频率误差模型。

这里建立质心算法像素频率误差模型：

A_x(sin(2πx_p+2πB_x)+A_y(sin(2πy_p+2πB_y)-sin(2πB_x)-sin(2πB_y))  (4)

这里Ax，Ay为能量误差系数，Bx，By为相位误差系数，C为水平误差系数，(x_p，y_p)为一个像素内的位置坐标。在整个靶面上，Ax，Ay为常数，Bx，By为2阶变化。对上述系数在整个靶面上进行建模为：

A_x＝A_y＝T

$B_x=\left(\frac{x}{N}\right)b_1-{\left(\frac{x}{N}\right)}^2{b}_2---\left(5\right)$

$B_y={\left(\frac{y}{N}\right)}^2{b}_1-{\left(\frac{y}{N}\right)}^2{b}_2$

这里，T为常数，N为靶面行或者列像素数，(x，y)为靶面像素整数坐标，b1，b2为相位偏差系数。

根据公式(4)我们可以得到像素频率偏差函数：

f(T，Bx，By)＝

T(sin(2πx_p+2πB_x)+T(sin(2πy_p+2πB_y)-sin(2πB_x)-sin(2πB_y))  (6)

$\mathrm{\&Delta;}d=\left(\begin{array}{ccc}\frac{\&PartialD;f}{\&PartialD;T}& \frac{\&PartialD;f}{\&PartialD;\mathrm{Bx}}& \frac{\&PartialD;f}{\&PartialD;\mathrm{By}}\end{array}\right)\left(\begin{array}{c}\mathrm{\&Delta;T}\\ \mathrm{\&Delta;Bx}\\ \mathrm{\&Delta;By}\end{array}\right)---\left(7\right)$

这里，Δd为参数估计误差，采用最小二乘迭代估计方法，可以求得以上3个参数的估计值。然后，根据计算得到的各个像素Bx和By的估计值，道理来计算整个靶面的相位偏差系数b1和b2。根据公式(5)有：

$\mathrm{\&Delta;}\mathrm{Bx}=\left[\begin{array}{cc}\frac{x}{N},& -{\left(\frac{x}{N}\right)}^2\end{array}\right]\left[\begin{array}{c}\mathrm{\&Delta;}{b}_1\\ \mathrm{\&Delta;}{b}_2\end{array}\right]---\left(8\right)$

$\mathrm{\&Delta;By}=\left[\begin{array}{cc}\frac{y}{N},& -{\left(\frac{y}{N}\right)}^2\end{array}\right]\left[\begin{array}{c}\mathrm{\&Delta;}{b}_1\\ \mathrm{\&Delta;}{b}_2\end{array}\right]$

然后，联立两方程采用最小二乘法就可以估计出b1和b2。

数据采集。

首先，考虑在一个像素内的采集，如图4所示，采集点阵为9×9点阵。根据9×9的点阵数据计算位置偏差，然后来提取该像素的模型参数。

这里采用的星敏感器靶面为1024×1204像素阵列，同样选择整个靶面的像素进行9×9采样，x方向采样位置为：

[1 129 257 385 513 641 769 897 1024]

同理y方向采样间距也为128。于是得到9×9的靶面采集像素位置，根据这些位置计算的b1和b2。

仿真和误差分析。

仿真的星敏感器基本参数为：

视场：12°×12°；

像素阵列：1024×1024；

像素尺寸：0.015m×0.015mm；

焦距：73.6059mm

假设像素频率内，星敏感器的随机误差为0.005，各个系数的设定值为：

T＝0.05；b1＝0.1667，b2＝0.1667；

以(1，1)像素为例，没有噪声的影响下，其Bx＝-0.2618，By＝-0.2618，下表表示了最小二乘法迭代求值过程：

	1	2	3	4	5
						T	0.048296	0.048296	0.049999	0.05	0.05
Bx	0	-0.26795	-0.26158	-0.2618	-0.2618
						By	0	-0.26795	-0.26158	-0.2618	-0.2618

如果添加了0.001像素高斯噪声，那么得到的计算结果如下：

	1	2	3	4	5
						T	0.048329	0.048329	0.05005	0.050051	0.050051
Bx	0	-0.27316	-0.26627	-0.26652	-0.26652
						By	0	-0.2654	-0.25936	-0.25956	-0.25956

将上述参数带入实际工作过程，验证所得残留校准误差均值为：0.0057像素，说明已经取得了很好的校准效果。

然后根据Bx和By值来估计参数b1和b2。在没有噪声的影响下，估计值为：b1＝1.0472，b2＝1.0472。假设噪声为0.01方差的高斯噪声，则估计结果为：b1＝1.0469，b2＝1.0522。

通过结果可以看出，经过该像素频率级的高精度星敏感器校准，则星敏感器的质心求取精度可以达到1/100～1/200水平，大大高于未经像素频率级校准的星敏感器。

Claims (1)

1.一种基于像素频率的星敏感器校准方法，使用由气垫平台、单星星光模拟器、星敏感器、2维轴向转台以及数据处理计算机组成的校准系统，其特征在于，

1.1、星敏感器质心算法；

1.1.1、建立点扩散函数模型，星像在图像传感器靶面上的分布称为点扩散函数，通过计算弥散圆斑的质心来获得亚像素精度的星点位置坐标，点扩散函数模型见式(1)：

$I(x,y)=\left\{\begin{array}{cc}{I}_0\exp (-(x^2+y^2))& x^2+y^2\&GreaterEqual;r^2\\ I_0& x^2+y^2<r^2\end{array}\right.---\left(1\right)$

式中，I₀为中心信号输出，x，y为像素坐标，r为截断半径，该模型是点扩散函数在某个高度截断后的形状；

1.1.2、建立质心算法公式为：

$x_c=\frac{{\&Integral;\&Integral;}_4\mathrm{xI}(x,y)\mathrm{dxdy}}{{\&Integral;\&Integral;}_{4}I(x,y)\mathrm{dxdy}},y_c=\frac{{\&Integral;\&Integral;}_4\mathrm{yI}(x,y)\mathrm{dxdy}}{{\&Integral;\&Integral;}_{4}I(x,y)\mathrm{dxdy}}---\left(2\right)$

式中，x_c，y_c为星点辐射中心，A表示大于阈值的区域，x，y为靶面像素坐标，I(x，y)为辐射分布函数，经过靶面像素的离散化后，质心计算公式为：

${\tilde{x}}_c=\frac{{\mathrm{\&Sigma;}}_{k=1}^n{x}_k{I}_k}{{\mathrm{\&Sigma;}}_{k=1}^n{I}_k},{\tilde{y}}_c=\frac{{\mathrm{\&Sigma;}}_{k=1}^n{y}_k{I}_k}{{\mathrm{\&Sigma;}}_{k=1}^n{I}_k}---\left(3\right)$

式中，

为离散化后的星点中心，n为阈值区域内像素个数，k为像素序号，x_k，y_k第k个像素的坐标，I_k为第k个像素内积累的光电子对应亮度输出；

1.2、建立质心算法像素频率误差模型；

1.2.1、质心算法像素频率误差模型见式(4)：

A_x(sin(2πx_p+2πB_x)+Ay(sin(2πy_p+2πB_y)-sin(2πB_x)-sin(2πB_y))    (4)

式中，Ax、Ay为能量误差系数，Bx，By为相位误差系数，(x_p，y_p)为一个像素内的位置坐标；在整个靶面上，Ax，Ay为常数，Bx，By为2阶变化；

1.2.2、对上述系数在整个靶面上进行建模为：

A_x＝A_y＝T

$B_x=\left(\frac{x}{N}\right)b_1-{\left(\frac{x}{N}\right)}^2{b}_2---\left(5\right)$

$B_y={\left(\frac{y}{N}\right)}^2{b}_1-{\left(\frac{y}{N}\right)}^2{b}_2$

式中，T为常数，N为靶面行或者列像素数，(x，y)为靶面像素整数坐标，b1，b2为相位偏差系数；

1.2.3、根据公式(4)得到像素频率偏差函数：

f(T，Bx，By)＝

T(sin(2πx_p+2πB_x)+T(sin(2πy_p+2πB_y)-sin(2πB_x)-sin(2πB_y))  (6)

$\mathrm{\&Delta;d}=\left(\begin{array}{ccc}\frac{\&PartialD;f}{\&PartialD;T}& \frac{\&PartialD;f}{\&PartialD;\mathrm{Bx}}& \frac{\&PartialD;f}{\&PartialD;\mathrm{By}}\end{array}\right)\left(\begin{array}{c}\mathrm{\&Delta;T}\\ \mathrm{\&Delta;Bx}\\ \mathrm{\&Delta;By}\end{array}\right)---\left(7\right)$

式中，Δd为参数估计误差，ΔT为能量误差，ΔBx、ΔBy分别为x和y方向相位偏差，采用最小二乘迭代估计方法，求得以上3个参数的估计值；然后计算得到各个像素Bx和By的估计值，并计算整个靶面的相位偏差系数b1和b2；根据公式(5)有：

$\mathrm{\&Delta;Bx}=\left[\begin{array}{cc}\frac{x}{N},& -{\left(\frac{x}{N}\right)}^2\end{array}\right]\left[\begin{array}{c}\mathrm{\&Delta;}{b}_1\\ \mathrm{\&Delta;}{b}_2\end{array}\right]---\left(8\right)$

$\mathrm{\&Delta;By}=\left[\begin{array}{cc}\frac{y}{N},& -{\left(\frac{y}{N}\right)}^2\end{array}\right]\left[\begin{array}{c}\mathrm{\&Delta;}{b}_1\\ \mathrm{\&Delta;}{b}_2\end{array}\right]$

然后，联立两方程采用最小二乘法估计出b1和b2；

1.3、数据采集；

1.3.1、在一个像素内的采集，采集点阵为9×9点阵；根据9×9的点阵数据计算位置偏差，然后来提取该像素的模型参数；

1.3.2、选择整个靶面的像素进行9×9采样，所采用的星敏感器靶面为1024×1204像素阵列，x方向采样位置为：1、129、257、385、513、641、769、897、1024，采样间距为128；在y方向采样间距也为128，得到9×9的靶面采集像素位置；将这些位置上的测量值和估计值之间的偏差代入公式(8)中计算得到b1和b2。

Images