第三章 太阳矢量确定：从理论到工程实现

各位同学，今天我们来聊聊卫星太阳翼对日定向中最核心的一环——怎么确定太阳在哪儿。

说白了，你要让太阳翼对准太阳，首先得知道太阳在哪个方向。这个看似简单的问题，在航天工程里其实有不少门道。我当年做第一个型号的时候，就因为在太阳矢量计算上栽了个跟头，导致整星在轨测试多花了两周时间。嗯，从那以后我对这块就特别上心。

3.1 太阳在惯性系下的位置计算

先说说最基础的问题：怎么知道太阳在惯性坐标系里的位置？

你可能会想，太阳不是基本不动吗？其实不是。地球绕着太阳转，卫星绕着地球转，这个相对运动必须算清楚。

工程上常用的方法是查表加插值。NASA 提供的 DE 历表（比如 DE405）精度很高，但星上计算机资源有限，不可能实时算那么复杂的模型。我个人习惯的做法是：

• 轨道预报阶段：用高精度模型算好太阳位置，做成查找表
• 在轨运行阶段：用简化解析公式实时计算

这里给出一套工程上常用的简化算法：

/* 太阳位置简化计算（基于J2000.0惯性系） */
/* 输入：MJD（简化儒略日） */
/* 输出：太阳单位方向向量 */

void sun_vector_calc(double mjd, double sun_vec[3])
{
    double T = (mjd - 51544.5) / 36525.0;  // 从J2000.0起算的儒略世纪数
    
    double M = 357.5291 + 35999.0503 * T;   // 太阳平近点角（度）
    double lambda = 280.4665 + 36000.7698 * T; // 太阳平黄经
    
    /* 中心差修正 */
    double C = (1.9146 - 0.0048 * T) * sin(M * DEG2RAD) 
             + 0.0200 * sin(2 * M * DEG2RAD);
    
    double true_lambda = lambda + C;  // 太阳真黄经
    
    /* 计算黄道坐标系下的位置 */
    double r = 1.000001 - 0.01671 * cos(M * DEG2RAD) 
             - 0.00014 * cos(2 * M * DEG2RAD);  // 日地距离（AU）
    
    double x = r * cos(true_lambda * DEG2RAD);
    double y = r * sin(true_lambda * DEG2RAD);
    
    /* 转换到J2000.0惯性系（考虑黄赤交角） */
    double eps = 23.4393 - 0.0130 * T;  // 黄赤交角
    
    sun_vec[0] = x;
    sun_vec[1] = y * cos(eps * DEG2RAD);
    sun_vec[2] = y * sin(eps * DEG2RAD);
    
    /* 归一化 */
    double norm = sqrt(sun_vec[0]*sun_vec[0] + 
                      sun_vec[1]*sun_vec[1] + 
                      sun_vec[2]*sun_vec[2]);
    sun_vec[0] /= norm;
    sun_vec[1] /= norm;
    sun_vec[2] /= norm;
}

3.2 太阳敏感器测量模型

理论计算是一回事，实际测量是另一回事。卫星上装什么传感器来测太阳方向？

常用的太阳敏感器分两类：

• 模拟式太阳敏感器：输出太阳角的正弦/余弦值，精度一般，但简单可靠
• 数字式太阳敏感器：用APS/CCD成像，直接输出太阳中心在像平面上的坐标，精度高

我做过一个项目，用的是数字式太阳敏感器。它的测量模型长这样：

/* 数字太阳敏感器测量模型 */
/* 输入：太阳在敏感器坐标系下的方向向量 s_sensor[3] */
/* 输出：像平面坐标 (u, v) 像素 */

void sun_sensor_model(double s_sensor[3], double *u, double *v)
{
    double f = 12.5e-3;  // 焦距（米）
    double pixel_size = 5.5e-6;  // 像元尺寸（米）
    
    /* 针孔相机模型 */
    double u_m = f * s_sensor[0] / s_sensor[2];  // 米
    double v_m = f * s_sensor[1] / s_sensor[2];
    
    /* 转换到像素坐标 */
    *u = u_m / pixel_size + 512;  // 假设1024x1024像元，中心在(512,512)
    *v = v_m / pixel_size + 512;
    
    /* 别忘了加噪声！ */
    /* 实际工程中还要考虑畸变、温度漂移等 */
}

3.3 矢量观测与滤波

有了理论计算的太阳位置，又有了太阳敏感器的实测数据，接下来要做的就是融合它们。

你想想看，理论计算值很平滑但可能有系统误差，实测值有噪声但能反映真实情况。怎么结合？

常用的方法有两种：

方法	优点	缺点	适用场景
最小二乘拟合	计算简单，适合资源受限	不能处理动态变化	稳态对日定向
扩展卡尔曼滤波（EKF）	能处理非线性，精度高	计算量大，调参麻烦	高精度指向任务
无迹卡尔曼滤波（UKF）	比EKF更稳定	计算量更大	强非线性场景

我个人比较推荐用EKF。为什么？因为它在精度和计算量之间取得了很好的平衡。下面是我在项目中用过的EKF实现框架：

/* 太阳矢量EKF滤波器（简化版） */
/* 状态量：太阳方向向量 (3维) */
/* 观测量：太阳敏感器输出 (2维：u, v) */

typedef struct {
    double x[3];      // 状态：太阳方向向量
    double P[3][3];   // 状态协方差矩阵
    double Q[3][3];   // 过程噪声
    double R[2][2];   // 测量噪声
} sun_ekf_t;

void sun_ekf_predict(sun_ekf_t *ekf, double dt)
{
    /* 预测步：太阳方向变化很慢，用单位矩阵 */
    /* 实际上要考虑卫星轨道运动引起的太阳方向变化 */
    double F[3][3] = {{1,0,0},{0,1,0},{0,0,1}};
    
    /* 状态预测：x = F * x */
    /* 协方差预测：P = F*P*F' + Q */
    /* 具体实现略... */
}

void sun_ekf_update(sun_ekf_t *ekf, double u_meas, double v_meas)
{
    /* 更新步：用太阳敏感器测量值修正 */
    /* 计算卡尔曼增益 */
    /* 更新状态和协方差 */
    /* 具体实现略... */
}

3.4 工程实现中的坑

最后说说我踩过的几个坑，希望对你有帮助：

• 太阳敏感器安装偏差：装上去之后一定要做地面标定，否则0.1度的安装误差就能让太阳翼效率下降不少
• 温度影响：太阳敏感器在轨温度变化很大，我见过因为温度导致测量值漂移0.5度的情况。解决办法是加温度补偿模型
• 地影区处理：卫星进入地影后，太阳敏感器没信号了。这时候要用理论计算值保持，等出影后再切换回测量值。切换时要平滑过渡，不然会有冲击

嗯，关于太阳矢量确定的内容就讲到这里。这部分虽然基础，但直接影响整个对日定向系统的性能。希望你能把这些方法真正用到自己的项目中去。