J4 Analytical Propagator
J4 Analytical Propagator
Note: This document is AI-translated.
Overview
The J4 analytical propagator provides orbit propagation considering Earth's J2 and J4 oblateness terms, offering higher precision orbit prediction than the J2 propagator. This method maintains high computational efficiency while accounting for more detailed Earth non-spherical gravitational field effects.
Core Concepts
J4 Perturbation Term
The J4 term is the second-order term in Earth's gravitational field spherical harmonic expansion, representing more detailed Earth shape characteristics. J4 perturbation further refines:
- Orbital plane precession: More accurate orbital plane precession
- Perigee precession: More precise perigee rotation
- Higher-order effects: Additional perturbation effects beyond J2
Propagation Principle
The J4 analytical propagator is based on long-term variation rates of mean orbital elements considering both J2 and J4 terms:
- Calculate long-term variation rates of orbital elements caused by J2 and J4 perturbations
- Perform linear propagation of mean orbital elements
- Convert mean orbital elements to instantaneous orbital elements
- Calculate corresponding position and velocity vectors
Main Interfaces
Class Interface
class J4Analytical : public J2J4AnalyticalConstructors:
J4Analytical(const ModOrbElem& modOrbElem, const TimePoint& epoch, double gm, double j2, double j4, double re);
J4Analytical(const Vector3d& pos, const Vector3d& vel, const TimePoint& epoch, double gm, double j2, double j4, double re);Function Interface
AST_CORE_API errc_t aJ4AnalyticalProp(double duration, double gm, double j2, double j4, double re, Vector3d& r, Vector3d& v);
AST_CORE_API errc_t aJ4AnalyticalProp(double duration, double gm, double j2, double j4, double re, ModOrbElem& modOrbElem);Parameter Description:
duration: Propagation time interval (seconds)gm: Gravitational constant (m³/s²)j2: J2 coefficient (dimensionless)j4: J4 coefficient (dimensionless)re: Body radius (meters)r: Input/output position vector (meters)v: Input/output velocity vector (m/s)modOrbElem: Input/output modified orbital elements
Usage Examples
1. Position-Velocity Propagation
#include "AstCore/J4Analytical.hpp"
#include "AstCore/TimePoint.hpp"
#include "AstMath/Vector.hpp"
#include "AstUtil/Literals.hpp"
#include "AstUtil/Constants.h"
#include <iostream>
#include <iomanip>
AST_USING_NAMESPACE
using namespace _AST literals;
int main()
{
// 地球参数
const double gm_earth = kEarthGrav;
const double j2_earth = 1.08262668e-3; // 地球J2系数
const double j4_earth = -1.620e-6; // 地球J4系数
const double re_earth = kEarthRadius; // 地球半径
// 初始轨道状态 (近地轨道)
Vector3d position{7000_km, 0.0, 0.0};
Vector3d velocity{0.0, 7.5_km_s, 1.0_km_s};
std::cout << std::fixed << std::setprecision(3);
std::cout << "J4解析传播器 - 位置速度传播示例" << std::endl;
std::cout << "==================================" << std::endl;
std::cout << "初始轨道状态:" << std::endl;
std::cout << "位置: (" << position[0]/1000.0 << ", "
<< position[1]/1000.0 << ", "
<< position[2]/1000.0 << ") km" << std::endl;
std::cout << "速度: (" << velocity[0]/1000.0 << ", "
<< velocity[1]/1000.0 << ", "
<< velocity[2]/1000.0 << ") km/s" << std::endl;
std::cout << std::endl;
// 传播时间间隔 (3天)
double duration = 3 * 24 * 3600.0; // 秒
// 执行J4解析传播
errc_t result = aJ4AnalyticalProp(duration, gm_earth, j2_earth, j4_earth, re_earth, position, velocity);
if (result == eNoError) {
std::cout << "传播后轨道状态 (" << duration/(24*3600.0) << " 天后):" << std::endl;
std::cout << "位置: (" << position[0]/1000.0 << ", "
<< position[1]/1000.0 << ", "
<< position[2]/1000.0 << ") km" << std::endl;
std::cout << "速度: (" << velocity[0]/1000.0 << ", "
<< velocity[1]/1000.0 << ", "
<< velocity[2]/1000.0 << ") km/s" << std::endl;
// 计算轨道倾角变化
double initial_inclination = std::atan2(1.0, 7.5) * kRadToDeg; // 近似计算
double final_inclination = std::atan2(velocity[2], std::sqrt(velocity[0]*velocity[0] + velocity[1]*velocity[1])) * kRadToDeg;
std::cout << std::endl;
std::cout << "轨道参数:" << std::endl;
std::cout << "初始倾角: " << initial_inclination << " 度" << std::endl;
std::cout << "最终倾角: " << final_inclination << " 度" << std::endl;
std::cout << "倾角变化: " << (final_inclination - initial_inclination) << " 度" << std::endl;
} else {
std::cout << "传播失败,错误码: " << result << std::endl;
}
return 0;
}2. Orbital Elements Propagation
#include "AstCore/J4Analytical.hpp"
#include "AstCore/TimePoint.hpp"
#include "AstCore/OrbitElement.hpp"
#include "AstMath/Vector.hpp"
#include "AstUtil/Literals.hpp"
#include "AstUtil/Constants.h"
#include <iostream>
#include <iomanip>
#include <cmath>
AST_USING_NAMESPACE
using namespace _AST literals;
int main()
{
// 地球参数
const double gm_earth = kEarthGrav;
const double j2_earth = 1.08262668e-3;
const double j4_earth = -1.620e-6;
const double re_earth = kEarthRadius;
// 创建改进轨道根数 (大椭圆轨道示例)
ModOrbElem orb_elem;
orb_elem.rp_ = 26500_km; // 近拱点半径
orb_elem.e_ = 0.73; // 偏心率
orb_elem.i_ = 28.5 * kDegToRad; // 倾角
orb_elem.raan_ = 0.0;
orb_elem.argper_ = 0.0;
orb_elem.trueA_ = 0.0;
std::cout << std::fixed << std::setprecision(6);
std::cout << "J4解析传播器 - 轨道根数传播示例" << std::endl;
std::cout << "==================================" << std::endl;
std::cout << "初始轨道根数 (大椭圆轨道):" << std::endl;
std::cout << "半长轴: " << orb_elem.getA() / 1000.0 << " km" << std::endl;
std::cout << "偏心率: " << orb_elem.e_ << std::endl;
std::cout << "倾角: " << orb_elem.i_ * kRadToDeg << " 度" << std::endl;
std::cout << "升交点赤经: " << orb_elem.raan_ * kRadToDeg << " 度" << std::endl;
std::cout << "近地点幅角: " << orb_elem.argper_ * kRadToDeg << " 度" << std::endl;
std::cout << std::endl;
// 传播多个时间点,观察J4摄动效应
double duration_day = 7.0; // 7天
double duration_sec = duration_day * 24 * 3600.0;
std::cout << "轨道根数随时间变化 (" << duration_day << " 天):" << std::endl;
std::cout << "时间(天)\t升交点赤经(度)\t近地点幅角(度)\t偏心率" << std::endl;
std::cout << "----------------------------------------------------------------" << std::endl;
for (int i = 0; i <= 7; i++) {
double time = i * duration_sec / 7.0;
// 复制初始轨道根数
ModOrbElem current_orb = orb_elem;
// 执行J4解析传播
errc_t result = aJ4AnalyticalProp(time, gm_earth, j2_earth, j4_earth, re_earth, current_orb);
if (result == eNoError) {
std::cout << std::setw(6) << time/(24*3600.0) << "\t"
<< std::setw(12) << current_orb.raan_ * kRadToDeg << "\t"
<< std::setw(12) << current_orb.argper_ * kRadToDeg << "\t"
<< std::setw(10) << current_orb.e_ << std::endl;
} else {
std::cout << "传播失败,错误码: " << result << std::endl;
break;
}
}
// 计算J4对近地点进动的修正
double semi_major_axis = orb_elem.getA();
double eccentricity = orb_elem.e_;
double inclination = orb_elem.i_;
// J2引起的近地点进动率
double n = std::sqrt(gm_earth / (semi_major_axis * semi_major_axis * semi_major_axis));
double p = semi_major_axis * (1 - eccentricity * eccentricity);
double aop_rate_j2 = 1.5 * j2_earth * n * std::pow(re_earth / p, 2) * (2 - 2.5 * std::sin(inclination) * std::sin(inclination));
// J4引起的近地点进动率修正
double aop_rate_j4 = 0.0; // 简化计算,实际需要复杂公式
std::cout << std::endl;
std::cout << "摄动理论值:" << std::endl;
std::cout << "J2近地点进动率: " << aop_rate_j2 * kRadToDeg * 86400.0 << " 度/天" << std::endl;
std::cout << "J4提供更高精度修正" << std::endl;
return 0;
}3. Class Interface Usage
#include "AstCore/J4Analytical.hpp"
#include "AstCore/TimePoint.hpp"
#include "AstCore/OrbitElement.hpp"
#include "AstMath/Vector.hpp"
#include "AstUtil/Literals.hpp"
#include "AstUtil/Constants.h"
#include <iostream>
#include <iomanip>
AST_USING_NAMESPACE
using namespace _AST literals;
int main()
{
std::cout << std::fixed << std::setprecision(6);
std::cout << "J4解析传播器 - 类接口使用示例" << std::endl;
std::cout << "================================" << std::endl;
// 地球参数
const double gm_earth = 3.986004418e14; // 地球引力常数 (m³/s²)
const double j2_earth = 1.08262668e-3; // 地球J2系数
const double j4_earth = -1.620e-6; // 地球J4系数
const double re_earth = 6378137.0; // 地球半径 (m)
// 示例1: 使用轨道根数构造J4传播器
std::cout << "示例1: 使用轨道根数构造J4传播器" << std::endl;
std::cout << "--------------------------------" << std::endl;
TimePoint epoch1 = TimePoint::FromUTC(2026, 1, 1, 0, 0, 0);
// 创建改进轨道根数 (太阳同步轨道)
ModOrbElem modOrbElem1;
modOrbElem1.rp_ = 7178.137_km; // 近拱点半径 (km)
modOrbElem1.e_ = 0.001; // 偏心率
modOrbElem1.i_ = 98.0_deg; // 轨道倾角
modOrbElem1.argper_ = 0.0_deg; // 近拱点角
modOrbElem1.raan_ = 0.0_deg; // 升交点赤经
modOrbElem1.trueA_ = 0.0_deg; // 真近点角
// 使用轨道根数构造J4传播器
J4Analytical j4_prop1(modOrbElem1, epoch1, gm_earth, j2_earth, j4_earth, re_earth);
std::cout << "轨道根数构造成功" << std::endl;
std::cout << "半长轴: " << modOrbElem1.getA()/1000.0 << " km" << std::endl;
std::cout << "偏心率: " << modOrbElem1.e_ << std::endl;
std::cout << "轨道倾角: " << modOrbElem1.i_/kDegToRad << " 度" << std::endl;
std::cout << std::endl;
// 示例2: 使用位置速度构造J4传播器
std::cout << "示例2: 使用位置速度构造J4传播器" << std::endl;
std::cout << "--------------------------------" << std::endl;
TimePoint epoch2 = TimePoint::FromUTC(2026, 1, 1, 0, 0, 0);
// 初始位置速度 (近圆轨道)
Vector3d position{7000_km, 0.0, 0.0};
Vector3d velocity{0.0, 7.5_km_s, 0.0};
// 使用位置速度构造J4传播器
J4Analytical j4_prop2(position, velocity, epoch2, gm_earth, j2_earth, j4_earth, re_earth);
std::cout << "位置速度构造成功" << std::endl;
std::cout << "位置: (" << position[0]/1000.0 << ", "
<< position[1]/1000.0 << ", "
<< position[2]/1000.0 << ") km" << std::endl;
std::cout << "速度: (" << velocity[0]/1000.0 << ", "
<< velocity[1]/1000.0 << ", "
<< velocity[2]/1000.0 << ") km/s" << std::endl;
std::cout << std::endl;
// 示例3: 完整传播流程演示
std::cout << "示例3: 完整传播流程演示" << std::endl;
std::cout << "------------------------" << std::endl;
// 创建时间点
TimePoint start_time = TimePoint::FromUTC(2026, 1, 1, 0, 0, 0);
TimePoint end_time = start_time + 24 * 3600.0; // 1天后
// 初始轨道状态 (大椭圆轨道)
Vector3d pos{26560_km, 0.0, 0.0}; // 约20000km高度
Vector3d vel{0.0, 3.87_km_s, 0.0}; // 约12小时周期
std::cout << "传播前状态:" << std::endl;
std::cout << "起始时间: 2026-01-01 00:00:00 UTC" << std::endl;
std::cout << "位置: (" << pos[0]/1000.0 << ", "
<< pos[1]/1000.0 << ", "
<< pos[2]/1000.0 << ") km" << std::endl;
std::cout << "速度: (" << vel[0]/1000.0 << ", "
<< vel[1]/1000.0 << ", "
<< vel[2]/1000.0 << ") km/s" << std::endl;
// 使用函数接口进行传播
errc_t result = aJ4AnalyticalProp(24 * 3600.0, gm_earth, j2_earth, j4_earth, re_earth, pos, vel);
if (result == eNoError) {
std::cout << std::endl;
std::cout << "传播后状态:" << std::endl;
std::cout << "结束时间: 2026-01-02 00:00:00 UTC" << std::endl;
std::cout << "位置: (" << pos[0]/1000.0 << ", "
<< pos[1]/1000.0 << ", "
<< pos[2]/1000.0 << ") km" << std::endl;
std::cout << "速度: (" << vel[0]/1000.0 << ", "
<< vel[1]/1000.0 << ", "
<< vel[2]/1000.0 << ") km/s" << std::endl;
// 计算轨道参数变化
double initial_radius = 26560.0; // km
double final_radius = pos.norm() / 1000.0;
double radius_change = final_radius - initial_radius;
std::cout << std::endl;
std::cout << "轨道参数变化:" << std::endl;
std::cout << "初始轨道半径: " << initial_radius << " km" << std::endl;
std::cout << "最终轨道半径: " << final_radius << " km" << std::endl;
std::cout << "轨道半径变化: " << radius_change << " km" << std::endl;
std::cout << "注: 变化主要由J2和J4摄动引起" << std::endl;
} else {
std::cout << "传播失败,错误码: " << result << std::endl;
}
std::cout << std::endl;
std::cout << "J4解析传播器类接口使用总结:" << std::endl;
std::cout << "1. 支持使用轨道根数或位置速度构造传播器" << std::endl;
std::cout << "2. 提供函数接口和类接口两种使用方式" << std::endl;
std::cout << "3. 考虑J2和J4摄动项,精度高于J2传播器" << std::endl;
std::cout << "4. 适用于需要较高精度的近地轨道传播" << std::endl;
return 0;
}Accuracy and Limitations
Accuracy Characteristics
- Calculation Accuracy: Considers J2 and J4 term long-term perturbations, better than J2 propagator
- Applicability: Suitable for near-Earth orbits requiring more detailed Earth shape effects
- Computational Efficiency: High, suitable for applications requiring higher accuracy than J2
Usage Limitations
- Only considers J2 and J4 term perturbations, ignores higher-order terms and other perturbation forces
- For high-precision orbit prediction requirements, recommend using numerical methods like HPOP
- Errors accumulate during long-term propagation due to model simplification
Earth J2/J4 Parameters
For Earth orbit calculations, recommended parameters:
- J2 coefficient: 1.08262668e-3
- J4 coefficient: -1.620e-6
- Earth radius: 6378137.0 meters
- Gravitational constant: 3.986004418e14 m³/s²
Accuracy Comparison
Compared to J2 propagator, J4 propagator provides:
- More accurate perigee precession rate calculation
- Better orbital plane precession prediction
- Improved accuracy for orbits with specific inclinations
Notes
- Input position and velocity vectors should be in the same inertial coordinate system
- J2/J4 coefficients and body radius should match the central body
- Propagation time interval should not be too long to avoid numerical error accumulation
- J4 effects are more significant for certain orbital inclinations and eccentricities
API Reference
///
/// @file J4Analytical.hpp
/// @brief ~
/// @details ~
/// @author axel
/// @date 2026-01-05
/// @copyright 版权所有 (C) 2026-present, ast项目.
///
/// ast项目(https://github.com/space-ast/ast)
/// 本项目基于 Apache 2.0 开源许可证分发。
/// 您可在遵守许可证条款的前提下使用、修改和分发本软件。
/// 许可证全文请见:
///
/// http://www.apache.org/licenses/LICENSE-2.0
///
/// 重要须知:
/// 软件按"现有状态"提供,无任何明示或暗示的担保条件。
/// 除非法律要求或书面同意,作者与贡献者不承担任何责任。
/// 使用本软件所产生的风险,需由您自行承担。
#pragma once
#include "AstGlobal.h"
#include "J2J4Analytical.hpp"
AST_NAMESPACE_BEGIN
/*!
@addtogroup Propagator
@{
*/
class ModOrbElem;
class AST_CORE_API J4Analytical: public J2J4Analytical
{
public:
/// @brief 构造函数
/// @details ~
/// @param modOrbElem 改进轨道根数
/// @param epoch 历元
/// @param gm 引力常数
/// @param j2 J2项
/// @param j4 J4项
/// @param re 天体半径
J4Analytical(const ModOrbElem& modOrbElem, const TimePoint& epoch, double gm, double j2, double j4, double re);
/// @brief 构造函数
/// @details ~
/// @param pos 位置向量
/// @param vel 速度向量
/// @param epoch 历元
/// @param gm 引力常数
/// @param j2 J2项
/// @param j4 J4项
/// @param re 天体半径
J4Analytical(const Vector3d& pos, const Vector3d& vel, const TimePoint& epoch, double gm, double j2, double j4, double re);
private:
// double j2_;
// double j4_;
// double re_;
};
/// @brief J4轨道预报
/// @details
/// @param duration 时间间隔
/// @param gm 引力常数
/// @param j2 J2项
/// @param j4 J4项
/// @param re 天体半径
/// @param r 位置向量
/// @param v 速度向量
/// @return 错误码
AST_CORE_API errc_t aJ4AnalyticalProp(double duration, double gm, double j2, double j4, double re, Vector3d& r, Vector3d& v);
/// @brief J4轨道预报
/// @details
/// @param duration 时间间隔
/// @param gm 引力常数
/// @param j2 J2项
/// @param j4 J4项
/// @param re 天体半径
/// @param modOrbElem 改进轨道根数
/// @return 错误码
AST_CORE_API errc_t aJ4AnalyticalProp(double duration, double gm, double j2, double j4, double re, ModOrbElem& modOrbElem);
/*! @} */
AST_NAMESPACE_END