1. 头文件和数据结构定义

1.1 头文件定义 (tle2kepler.h)

#ifndef TLE2KEPLER_H
#define TLE2KEPLER_H

#include <iostream>
#include <string>
#include <vector>
#include <cmath>
#include <sstream>
#include <iomanip>
#include <stdexcept>

// 数学常数定义
const double PI = 3.14159265358979323846;
const double TWO_PI = 2.0 * PI;
const double DEG_TO_RAD = PI / 180.0;
const double RAD_TO_DEG = 180.0 / PI;
const double MU = 398600.4418;  // 地球引力常数 (km^3/s^2)
const double J2 = 1.08262668e-3;  // 地球扁率J2项
const double EARTH_RADIUS = 6378.137;  // 地球赤道半径 (km)
const double OMEGA_EARTH = 7.2921151467e-5;  // 地球自转角速度 (rad/s)

// 轨道六根数结构体
struct KeplerElements {
    double semi_major_axis;       // 半长轴 a (km)
    double eccentricity;          // 偏心率 e
    double inclination;           // 轨道倾角 i (度)
    double raan;                  // 升交点赤经 Ω (度)
    double argument_of_perigee;   // 近地点幅角 ω (度)
    double true_anomaly;          // 真近点角 ν (度)
    double mean_anomaly;          // 平近点角 M (度)
    double epoch;                 // 历元 (儒略日)
    
    // 辅助参数
    double mean_motion;           // 平均运动 (rev/day)
    double bstar;                 // 阻力系数
    int revolution_number;        // 在轨圈数
    int satellite_catalog_number; // 卫星目录号
};

// TLE两行元素结构体
struct TLE {
    std::string line1;
    std::string line2;
    std::string name;  // 卫星名称
    
    // 解析后的数据
    int satellite_number;     // 卫星编号
    char classification;      // 分类等级
    int launch_year;         // 发射年份
    int launch_number;       // 发射序号
    std::string launch_piece; // 发射碎片编号
    int epoch_year;          // 历元年
    double epoch_day;        // 历元日
    double mean_motion_dot;  // 平均运动一阶导数
    double mean_motion_ddot; // 平均运动二阶导数
    double bstar_drag;       // BSTAR阻力系数
    int ephemeris_type;      // 星历类型
    int element_number;      // 元素编号
    double inclination;      // 轨道倾角 (度)
    double raan;             // 升交点赤经 (度)
    double eccentricity;     // 偏心率
    double argument_of_perigee;  // 近地点幅角 (度)
    double mean_anomaly;     // 平近点角 (度)
    double mean_motion;      // 平均运动 (rev/day)
    int revolution_number;   // 在轨圈数
};

// 主转换类
class TLE2Kepler {
private:
    // 内部转换函数
    static double calculateSemiMajorAxis(double mean_motion);
    static double meanToTrueAnomaly(double mean_anomaly, double eccentricity);
    static double eccentricToTrueAnomaly(double eccentric_anomaly, double eccentricity);
    static double solveKeplerEquation(double mean_anomaly, double eccentricity);
    
    // 解析辅助函数
    static double parseFloat(const std::string& str, int start, int length, 
                            double defaultValue = 0.0, double scale = 1.0);
    static int parseInt(const std::string& str, int start, int length, 
                       int defaultValue = 0);
    
public:
    // 解析TLE字符串
    static bool parseTLE(const std::string& tleLine1, 
                        const std::string& tleLine2, 
                        TLE& tleData);
    
    static bool parseTLE(const std::string& tleName,
                        const std::string& tleLine1, 
                        const std::string& tleLine2, 
                        TLE& tleData);
    
    // 从TLE计算开普勒根数
    static KeplerElements convertTLEtoKepler(const TLE& tleData);
    
    // 从TLE字符串直接计算
    static KeplerElements convertTLEtoKepler(const std::string& line1, 
                                            const std::string& line2);
    
    // 验证TLE格式
    static bool validateTLE(const std::string& line1, const std::string& line2);
    
    // 计算平均运动
    static double calculateMeanMotion(double semi_major_axis);
    
    // 轨道周期计算
    static double calculateOrbitalPeriod(double semi_major_axis);
    
    // 从平均运动计算半长轴
    static double meanMotionToSemiMajorAxis(double mean_motion);
    
    // 输出格式化
    static void printKeplerElements(const KeplerElements& kepler);
    static void printTLEData(const TLE& tle);
    
    // 坐标转换工具函数
    static void keplerToECI(const KeplerElements& kepler, double time,
                           double& x, double& y, double& z,
                           double& vx, double& vy, double& vz);
};

2. TLE解析实现

2.1 TLE解析实现 (tle2kepler.cpp)

#include "tle2kepler.h"

// 解析浮点数
double TLE2Kepler::parseFloat(const std::string& str, int start, int length, 
                              double defaultValue, double scale) {
    if (start < 0 || start + length > static_cast<int>(str.length())) {
        return defaultValue;
    }
    
    std::string sub = str.substr(start, length);
    
    // 去除前导空格
    size_t pos = sub.find_first_not_of(' ');
    if (pos == std::string::npos) {
        return defaultValue;
    }
    
    sub = sub.substr(pos);
    
    // 检查是否为空
    if (sub.empty()) {
        return defaultValue;
    }
    
    // 处理指数表示
    size_t e_pos = sub.find('-', 1);
    if (e_pos != std::string::npos && e_pos > 0 && 
        sub[e_pos-1] != 'E' && sub[e_pos-1] != 'e') {
        sub.insert(e_pos, "E");
    }
    
    try {
        double value = std::stod(sub) * scale;
        return value;
    } catch (const std::exception&) {
        return defaultValue;
    }
}

// 解析整数
int TLE2Kepler::parseInt(const std::string& str, int start, int length, 
                        int defaultValue) {
    if (start < 0 || start + length > static_cast<int>(str.length())) {
        return defaultValue;
    }
    
    std::string sub = str.substr(start, length);
    
    // 去除前导空格
    size_t pos = sub.find_first_not_of(' ');
    if (pos == std::string::npos) {
        return defaultValue;
    }
    
    sub = sub.substr(pos);
    
    try {
        return std::stoi(sub);
    } catch (const std::exception&) {
        return defaultValue;
    }
}

// 解析TLE
bool TLE2Kepler::parseTLE(const std::string& tleLine1, 
                         const std::string& tleLine2, 
                         TLE& tleData) {
    return parseTLE("", tleLine1, tleLine2, tleData);
}

bool TLE2Kepler::parseTLE(const std::string& tleName,
                         const std::string& tleLine1, 
                         const std::string& tleLine2, 
                         TLE& tleData) {
    // 验证行长度
    if (tleLine1.length() < 69 || tleLine2.length() < 69) {
        std::cerr << "Error: TLE lines are too short" << std::endl;
        return false;
    }
    
    tleData.name = tleName;
    tleData.line1 = tleLine1;
    tleData.line2 = tleLine2;
    
    // 解析第一行
    tleData.satellite_number = parseInt(tleLine1, 2, 5);
    tleData.classification = tleLine1[7];
    tleData.launch_year = parseInt(tleLine1, 9, 2);
    tleData.launch_number = parseInt(tleLine1, 11, 3);
    
    // 发射碎片编号
    std::string launch_piece_str = tleLine1.substr(14, 3);
    size_t pos = launch_piece_str.find_last_not_of(' ');
    tleData.launch_piece = (pos != std::string::npos) ? 
                           launch_piece_str.substr(0, pos + 1) : "";
    
    // 历元年份
    tleData.epoch_year = parseInt(tleLine1, 18, 2);
    tleData.epoch_day = parseFloat(tleLine1, 20, 12);
    
    // 平均运动一阶导数
    tleData.mean_motion_dot = parseFloat(tleLine1, 33, 10, 0.0, 2.0);
    
    // 平均运动二阶导数
    std::string mdd_str = tleLine1.substr(44, 6);
    double mdd_value = parseFloat(tleLine1, 44, 6, 0.0, 1.0);
    
    // 检查指数
    std::string exp_str = tleLine1.substr(50, 2);
    int exponent = parseInt(exp_str, 0, 2, 0);
    
    if (mdd_str[0] == '-') {
        tleData.mean_motion_ddot = -mdd_value * std::pow(10.0, exponent);
    } else {
        tleData.mean_motion_ddot = mdd_value * std::pow(10.0, exponent);
    }
    
    // BSTAR阻力系数
    std::string bstar_str = tleLine1.substr(53, 6);
    double bstar_value = parseFloat(tleLine1, 53, 6, 0.0, 1.0);
    
    std::string bstar_exp_str = tleLine1.substr(59, 2);
    int bstar_exponent = parseInt(bstar_exp_str, 0, 2, 0);
    
    if (bstar_str[0] == '-') {
        tleData.bstar_drag = -bstar_value * std::pow(10.0, bstar_exponent);
    } else {
        tleData.bstar_drag = bstar_value * std::pow(10.0, bstar_exponent);
    }
    
    tleData.ephemeris_type = parseInt(tleLine1, 62, 1);
    tleData.element_number = parseInt(tleLine1, 64, 4);
    
    // 解析第二行
    tleData.inclination = parseFloat(tleLine2, 8, 8);
    tleData.raan = parseFloat(tleLine2, 17, 8);
    
    // 偏心率（需要除以10000000）
    tleData.eccentricity = parseFloat(tleLine2, 26, 7, 0.0, 1e-7);
    
    tleData.argument_of_perigee = parseFloat(tleLine2, 34, 8);
    tleData.mean_anomaly = parseFloat(tleLine2, 43, 8);
    tleData.mean_motion = parseFloat(tleLine2, 52, 11);
    tleData.revolution_number = parseInt(tleLine2, 63, 5);
    
    return true;
}

// 验证TLE格式
bool TLE2Kepler::validateTLE(const std::string& line1, const std::string& line2) {
    if (line1.length() < 69 || line2.length() < 69) {
        return false;
    }
    
    if (line1[0] != '1' || line2[0] != '2') {
        return false;
    }
    
    // 验证校验和
    auto calculateChecksum = const std::string& line -> int {
        int sum = 0;
        for (size_t i = 0; i < 68; i++) {
            char c = line[i];
            if (c >= '0' && c <= '9') {
                sum += (c - '0');
            } else if (c == '-') {
                sum += 1;
            }
        }
        return sum % 10;
    };
    
    int checksum1 = parseInt(line1, 68, 1);
    int checksum2 = parseInt(line2, 68, 1);
    
    if (calculateChecksum(line1) != checksum1) {
        std::cerr << "Warning: Line 1 checksum mismatch" << std::endl;
    }
    
    if (calculateChecksum(line2) != checksum2) {
        std::cerr << "Warning: Line 2 checksum mismatch" << std::endl;
    }
    
    return true;
}

3. 开普勒方程求解

3.1 开普勒方程求解 (kepler_solver.cpp)

#include "tle2kepler.h"

// 求解开普勒方程 M = E - e*sin(E)
double TLE2Kepler::solveKeplerEquation(double mean_anomaly, double eccentricity) {
    // 转换为弧度
    double M = mean_anomaly * DEG_TO_RAD;
    double e = eccentricity;
    
    if (e < 0.0 || e >= 1.0) {
        return M;  // 无效偏心率，返回平近点角
    }
    
    // 初始猜测：E = M
    double E = M;
    
    // 牛顿-拉夫森迭代
    double tolerance = 1e-12;
    int max_iterations = 50;
    
    for (int i = 0; i < max_iterations; i++) {
        double f = E - e * sin(E) - M;
        double f_prime = 1.0 - e * cos(E);
        
        if (std::abs(f_prime) < 1e-12) {
            break;  // 避免除零
        }
        
        double delta = f / f_prime;
        E -= delta;
        
        if (std::abs(delta) < tolerance) {
            break;
        }
    }
    
    return E;  // 返回偏近点角（弧度）
}

// 从平近点角计算真近点角
double TLE2Kepler::meanToTrueAnomaly(double mean_anomaly, double eccentricity) {
    if (eccentricity < 0.0 || eccentricity >= 1.0) {
        return mean_anomaly;  // 无效偏心率
    }
    
    // 求解开普勒方程得到偏近点角
    double E = solveKeplerEquation(mean_anomaly, eccentricity);
    
    // 从偏近点角计算真近点角
    return eccentricToTrueAnomaly(E, eccentricity);
}

// 从偏近点角计算真近点角
double TLE2Kepler::eccentricToTrueAnomaly(double eccentric_anomaly, double eccentricity) {
    if (eccentricity < 0.0 || eccentricity >= 1.0) {
        return eccentric_anomaly * RAD_TO_DEG;
    }
    
    double e = eccentricity;
    double E = eccentric_anomaly;
    
    // 计算真近点角
    double sinE = sin(E);
    double cosE = cos(E);
    
    double sin_nu = sqrt(1.0 - e*e) * sinE / (1.0 - e * cosE);
    double cos_nu = (cosE - e) / (1.0 - e * cosE);
    
    double nu = atan2(sin_nu, cos_nu);
    
    // 确保角度在[0, 2π)范围内
    if (nu < 0.0) {
        nu += TWO_PI;
    }
    
    return nu * RAD_TO_DEG;  // 转换为度
}

4. 轨道参数计算

4.1 轨道参数计算 (orbital_calculations.cpp)

#include "tle2kepler.h"

// 从平均运动计算半长轴
double TLE2Kepler::calculateSemiMajorAxis(double mean_motion) {
    if (mean_motion <= 0.0) {
        return 0.0;
    }
    
    // 平均运动单位：转/天，需要转换为弧度/秒
    double n = mean_motion * TWO_PI / 86400.0;  // 弧度/秒
    
    // 开普勒第三定律: a^3 = μ / n^2
    double a = std::pow(MU / (n * n), 1.0/3.0);
    
    return a;
}

// 从平均运动计算半长轴（另一种方法）
double TLE2Kepler::meanMotionToSemiMajorAxis(double mean_motion) {
    if (mean_motion <= 0.0) {
        return 0.0;
    }
    
    // 更精确的考虑地球扁率的公式
    double n = mean_motion * TWO_PI / 86400.0;  // 转换为弧度/秒
    
    // 考虑J2项的一阶近似
    double a0 = std::pow(MU / (n * n), 1.0/3.0);
    
    // 迭代求解考虑J2影响的半长轴
    double a = a0;
    double tolerance = 1e-6;
    int max_iterations = 10;
    
    for (int i = 0; i < max_iterations; i++) {
        double p = a * (1.0 - 0.0);  // 简化，假设偏心率为0
        double n0 = sqrt(MU / (a * a * a));
        
        // J2项修正
        double delta_n = 1.5 * J2 * (EARTH_RADIUS * EARTH_RADIUS) * n0 * 
                        (1.0 - 0.0 * 0.0) / (p * p);  // 简化倾角
        
        double n_corrected = n0 + delta_n;
        double a_new = std::pow(MU / (n_corrected * n_corrected), 1.0/3.0);
        
        if (std::abs(a_new - a) < tolerance) {
            a = a_new;
            break;
        }
        
        a = a_new;
    }
    
    return a;
}

// 计算平均运动
double TLE2Kepler::calculateMeanMotion(double semi_major_axis) {
    if (semi_major_axis <= 0.0) {
        return 0.0;
    }
    
    double n = sqrt(MU / (semi_major_axis * semi_major_axis * semi_major_axis));
    
    // 转换为转/天
    return n * 86400.0 / TWO_PI;
}

// 计算轨道周期
double TLE2Kepler::calculateOrbitalPeriod(double semi_major_axis) {
    if (semi_major_axis <= 0.0) {
        return 0.0;
    }
    
    double period = TWO_PI * sqrt(semi_major_axis * semi_major_axis * semi_major_axis / MU);
    
    return period;  // 秒
}

5. 主转换函数

5.1 转换函数实现 (conversion.cpp)

#include "tle2kepler.h"

// 主转换函数：从TLE到开普勒根数
KeplerElements TLE2Kepler::convertTLEtoKepler(const TLE& tleData) {
    KeplerElements kepler;
    
    // 基本参数
    kepler.satellite_catalog_number = tleData.satellite_number;
    kepler.mean_motion = tleData.mean_motion;
    kepler.inclination = tleData.inclination;
    kepler.raan = tleData.raan;
    kepler.eccentricity = tleData.eccentricity;
    kepler.argument_of_perigee = tleData.argument_of_perigee;
    kepler.mean_anomaly = tleData.mean_anomaly;
    kepler.revolution_number = tleData.revolution_number;
    kepler.bstar = tleData.bstar_drag;
    
    // 计算历元（儒略日简化计算）
    int year = tleData.epoch_year;
    if (year < 57) {  // 1957-2056
        year += 2000;
    } else {
        year += 1900;
    }
    
    // 简化计算：年份 + 年积日
    kepler.epoch = static_cast<double>(year) + tleData.epoch_day;
    
    // 计算半长轴
    kepler.semi_major_axis = calculateSemiMajorAxis(tleData.mean_motion);
    
    // 计算真近点角
    kepler.true_anomaly = meanToTrueAnomaly(tleData.mean_anomaly, tleData.eccentricity);
    
    return kepler;
}

// 直接从TLE字符串转换
KeplerElements TLE2Kepler::convertTLEtoKepler(const std::string& line1, 
                                             const std::string& line2) {
    TLE tleData;
    
    if (!parseTLE(line1, line2, tleData)) {
        throw std::runtime_error("Failed to parse TLE");
    }
    
    return convertTLEtoKepler(tleData);
}

6. 坐标转换和输出

6.1 坐标转换 (coordinate_conversion.cpp)

#include "tle2kepler.h"

// 开普勒根数转地心惯性坐标系（简化版本）
void TLE2Kepler::keplerToECI(const KeplerElements& kepler, double time,
                            double& x, double& y, double& z,
                            double& vx, double& vy, double& vz) {
    // 提取开普勒根数
    double a = kepler.semi_major_axis;
    double e = kepler.eccentricity;
    double i = kepler.inclination * DEG_TO_RAD;
    double Omega = kepler.raan * DEG_TO_RAD;
    double omega = kepler.argument_of_perigee * DEG_TO_RAD;
    double nu = kepler.true_anomaly * DEG_TO_RAD;
    
    // 计算参数
    double p = a * (1.0 - e * e);  // 半通径
    
    if (p <= 0.0) {
        x = y = z = vx = vy = vz = 0.0;
        return;
    }
    
    // 计算位置在轨道平面内的坐标
    double r = p / (1.0 + e * cos(nu));
    double xp = r * cos(nu);
    double yp = r * sin(nu);
    
    // 计算速度在轨道平面内的分量
    double mu = MU;
    double h = sqrt(mu * p);
    double vxp = -mu / h * sin(nu);
    double vyp = mu / h * (e + cos(nu));
    
    // 旋转矩阵
    double cos_omega = cos(omega);
    double sin_omega = sin(omega);
    double cos_Omega = cos(Omega);
    double sin_Omega = sin(Omega);
    double cos_i = cos(i);
    double sin_i = sin(i);
    
    // 位置坐标转换
    x = (cos_omega * cos_Omega - sin_omega * sin_Omega * cos_i) * xp +
        (-sin_omega * cos_Omega - cos_omega * sin_Omega * cos_i) * yp;
    
    y = (cos_omega * sin_Omega + sin_omega * cos_Omega * cos_i) * xp +
        (-sin_omega * sin_Omega + cos_omega * cos_Omega * cos_i) * yp;
    
    z = (sin_omega * sin_i) * xp + (cos_omega * sin_i) * yp;
    
    // 速度坐标转换
    vx = (cos_omega * cos_Omega - sin_omega * sin_Omega * cos_i) * vxp +
         (-sin_omega * cos_Omega - cos_omega * sin_Omega * cos_i) * vyp;
    
    vy = (cos_omega * sin_Omega + sin_omega * cos_Omega * cos_i) * vxp +
         (-sin_omega * sin_Omega + cos_omega * cos_Omega * cos_i) * vyp;
    
    vz = (sin_omega * sin_i) * vxp + (cos_omega * sin_i) * vyp;
}

6.2 输出格式化 (output_formatter.cpp)

#include "tle2kepler.h"

// 输出开普勒根数
void TLE2Kepler::printKeplerElements(const KeplerElements& kepler) {
    std::cout << std::fixed << std::setprecision(6);
    std::cout << "\n========== Kepler Orbital Elements ==========\n";
    std::cout << "Satellite Catalog Number: " << kepler.satellite_catalog_number << "\n";
    std::cout << "Epoch: " << kepler.epoch << " (simplified JD)\n";
    std::cout << "Semi-major axis (a): " << kepler.semi_major_axis << " km\n";
    std::cout << "Eccentricity (e): " << kepler.eccentricity << "\n";
    std::cout << "Inclination (i): " << kepler.inclination << " deg\n";
    std::cout << "Right Ascension of Ascending Node (Ω): " << kepler.raan << " deg\n";
    std::cout << "Argument of Perigee (ω): " << kepler.argument_of_perigee << " deg\n";
    std::cout << "True Anomaly (ν): " << kepler.true_anomaly << " deg\n";
    std::cout << "Mean Anomaly (M): " << kepler.mean_anomaly << " deg\n";
    std::cout << "Mean Motion: " << kepler.mean_motion << " rev/day\n";
    std::cout << "BSTAR Drag Coefficient: " << kepler.bstar << "\n";
    std::cout << "Revolution Number: " << kepler.revolution_number << "\n";
    
    // 计算并显示轨道周期
    double period = calculateOrbitalPeriod(kepler.semi_major_axis);
    double period_min = period / 60.0;
    std::cout << "Orbital Period: " << period << " s (" << period_min << " min)\n";
    
    // 计算轨道高度
    double perigee = kepler.semi_major_axis * (1.0 - kepler.eccentricity) - EARTH_RADIUS;
    double apogee = kepler.semi_major_axis * (1.0 + kepler.eccentricity) - EARTH_RADIUS;
    std::cout << "Perigee Height: " << perigee << " km\n";
    std::cout << "Apogee Height: " << apogee << " km\n";
    
    std::cout << "==============================================\n";
}

// 输出TLE数据
void TLE2Kepler::printTLEData(const TLE& tle) {
    std::cout << std::fixed << std::setprecision(8);
    std::cout << "\n========== Parsed TLE Data ==========\n";
    if (!tle.name.empty()) {
        std::cout << "Satellite Name: " << tle.name << "\n";
    }
    std::cout << "Satellite Number: " << tle.satellite_number << "\n";
    std::cout << "Classification: " << tle.classification << "\n";
    std::cout << "Launch Year: " << tle.launch_year << "\n";
    std::cout << "Launch Number: " << tle.launch_number << "\n";
    std::cout << "Launch Piece: " << tle.launch_piece << "\n";
    std::cout << "Epoch Year: " << tle.epoch_year << "\n";
    std::cout << "Epoch Day: " << tle.epoch_day << "\n";
    std::cout << "Mean Motion (n): " << tle.mean_motion << " rev/day\n";
    std::cout << "Mean Motion First Derivative (n'): " << tle.mean_motion_dot << "\n";
    std::cout << "Mean Motion Second Derivative (n''): " << tle.mean_motion_ddot << "\n";
    std::cout << "BSTAR Drag: " << tle.bstar_drag << "\n";
    std::cout << "Inclination (i): " << tle.inclination << " deg\n";
    std::cout << "RAAN (Ω): " << tle.raan << " deg\n";
    std::cout << "Eccentricity (e): " << tle.eccentricity << "\n";
    std::cout << "Argument of Perigee (ω): " << tle.argument_of_perigee << " deg\n";
    std::cout << "Mean Anomaly (M): " << tle.mean_anomaly << " deg\n";
    std::cout << "Revolution Number: " << tle.revolution_number << "\n";
    std::cout << "======================================\n";
}

参考代码 tle转换为六根数的c++源代码 www.youwenfan.com/contentcsv/71159.html

7. 实用工具函数

7.1 儒略日计算 (julian_date.cpp)

#include "tle2kepler.h"

// 儒略日计算工具
class JulianDate {
public:
    // 计算儒略日
    static double toJulian(int year, int month, int day, 
                          int hour = 0, int minute = 0, 
                          double second = 0.0) {
        if (month <= 2) {
            year -= 1;
            month += 12;
        }
        
        int A = year / 100;
        int B = 2 - A + A / 4;
        
        double JD = static_cast<int>(365.25 * (year + 4716)) +
                   static_cast<int>(30.6001 * (month + 1)) +
                   day + B - 1524.5;
        
        // 添加时间
        double time = hour + minute / 60.0 + second / 3600.0;
        JD += time / 24.0;
        
        return JD;
    }
    
    // 从儒略日计算年月日
    static void fromJulian(double JD, int& year, int& month, int& day,
                          int& hour, int& minute, double& second) {
        JD += 0.5;
        int Z = static_cast<int>(JD);
        double F = JD - Z;
        
        int A = 0;
        if (Z < 2299161) {
            A = Z;
        } else {
            int alpha = static_cast<int>((Z - 1867216.25) / 36524.25);
            A = Z + 1 + alpha - alpha / 4;
        }
        
        int B = A + 1524;
        int C = static_cast<int>((B - 122.1) / 365.25);
        int D = static_cast<int>(365.25 * C);
        int E = static_cast<int>((B - D) / 30.6001);
        
        day = B - D - static_cast<int>(30.6001 * E);
        month = (E < 14) ? E - 1 : E - 13;
        year = (month > 2) ? C - 4716 : C - 4715;
        
        // 计算时间
        double time = F * 24.0;
        hour = static_cast<int>(time);
        time = (time - hour) * 60.0;
        minute = static_cast<int>(time);
        second = (time - minute) * 60.0;
    }
};

8. 示例和使用程序

8.1 主程序示例 (main.cpp)

#include "tle2kepler.h"
#include <fstream>
#include <vector>

// 示例TLE数据
const std::string EXAMPLE_TLE_NAME = "ISS (ZARYA)";
const std::string EXAMPLE_TLE_LINE1 = "1 25544U 98067A   22001.00000000  .00000000  00000-0  00000-0 0  9999";
const std::string EXAMPLE_TLE_LINE2 = "2 25544  51.6416  33.2520 0006926 359.9521  71.0405 15.49764999300001";

// 国际空间站TLE
const std::string ISS_TLE_NAME = "ISS (ZARYA)";
const std::string ISS_TLE_LINE1 = "1 25544U 98067A   22001.00000000  .00000000  00000-0  00000-0 0  9999";
const std::string ISS_TLE_LINE2 = "2 25544  51.6416  33.2520 0006926 359.9521  71.0405 15.49764999300001";

// GPS卫星TLE
const std::string GPS_TLE_NAME = "GPS BIIR-2  (PRN 13)";
const std::string GPS_TLE_LINE1 = "1 24876U 97035A   22001.00000000  .00000000  00000-0  00000-0 0  9999";
const std::string GPS_TLE_LINE2 = "2 24876  55.8093 266.9644 0190932 272.4905  85.8105  2.00562541123456";

// 测试函数
void testTLEConversion() {
    std::cout << "Testing TLE to Kepler Elements Conversion\n";
    std::cout << "=========================================\n";
    
    // 测试ISS
    std::cout << "\n1. International Space Station (ISS)\n";
    TLE issTLE;
    TLE2Kepler::parseTLE(ISS_TLE_NAME, ISS_TLE_LINE1, ISS_TLE_LINE2, issTLE);
    TLE2Kepler::printTLEData(issTLE);
    
    KeplerElements issKepler = TLE2Kepler::convertTLEtoKepler(issTLE);
    TLE2Kepler::printKeplerElements(issKepler);
    
    // 测试GPS卫星
    std::cout << "\n2. GPS Satellite\n";
    TLE gpsTLE;
    TLE2Kepler::parseTLE(GPS_TLE_NAME, GPS_TLE_LINE1, GPS_TLE_LINE2, gpsTLE);
    TLE2Kepler::printTLEData(gpsTLE);
    
    KeplerElements gpsKepler = TLE2Kepler::convertTLEtoKepler(gpsTLE);
    TLE2Kepler::printKeplerElements(gpsKepler);
    
    // 测试开普勒方程
    std::cout << "\n3. Testing Kepler Equation Solver\n";
    double test_mean_anomaly = 180.0;  // 度
    double test_eccentricity = 0.1;
    double true_anomaly = TLE2Kepler::meanToTrueAnomaly(test_mean_anomaly, test_eccentricity);
    std::cout << "Mean Anomaly: " << test_mean_anomaly << " deg\n";
    std::cout << "Eccentricity: " << test_eccentricity << "\n";
    std::cout << "True Anomaly: " << true_anomaly << " deg\n";
    
    // 测试半长轴计算
    std::cout << "\n4. Testing Semi-major Axis Calculation\n";
    double test_mean_motion = 15.0;  // rev/day
    double semi_major_axis = TLE2Kepler::calculateSemiMajorAxis(test_mean_motion);
    double orbital_period = TLE2Kepler::calculateOrbitalPeriod(semi_major_axis);
    std::cout << "Mean Motion: " << test_mean_motion << " rev/day\n";
    std::cout << "Semi-major Axis: " << semi_major_axis << " km\n";
    std::cout << "Orbital Period: " << orbital_period << " s (" 
              << orbital_period/60.0 << " min)\n";
}

// 从文件读取TLE
std::vector<TLE> readTLEFile(const std::string& filename) {
    std::vector<TLE> tle_list;
    std::ifstream file(filename);
    
    if (!file.is_open()) {
        std::cerr << "Error: Could not open file " << filename << std::endl;
        return tle_list;
    }
    
    std::string line;
    std::string name, line1, line2;
    
    while (std::getline(file, line)) {
        if (line.empty()) continue;
        
        if (line[0] == '1' && line.length() > 10) {
            line1 = line;
        } else if (line[0] == '2' && line.length() > 10) {
            line2 = line;
            
            if (!line1.empty()) {
                TLE tle;
                TLE2Kepler::parseTLE(name, line1, line2, tle);
                tle_list.push_back(tle);
                
                // 重置
                name.clear();
                line1.clear();
                line2.clear();
            }
        } else {
            // 这可能是卫星名称行
            name = line;
        }
    }
    
    file.close();
    return tle_list;
}

// 批量转换
void batchConversion(const std::string& filename) {
    std::vector<TLE> tle_list = readTLEFile(filename);
    
    std::cout << "\nBatch Conversion Results\n";
    std::cout << "=======================\n";
    std::cout << "Found " << tle_list.size() << " TLEs in file.\n\n";
    
    for (size_t i = 0; i < tle_list.size(); i++) {
        const TLE& tle = tle_list[i];
        
        std::cout << "\n[" << (i+1) << "] " << tle.name << "\n";
        
        try {
            KeplerElements kepler = TLE2Kepler::convertTLEtoKepler(tle);
            
            std::cout << "  Semi-major axis: " << kepler.semi_major_axis << " km\n";
            std::cout << "  Eccentricity: " << kepler.eccentricity << "\n";
            std::cout << "  Inclination: " << kepler.inclination << " deg\n";
            std::cout << "  Type: ";
            
            // 判断轨道类型
            double perigee = kepler.semi_major_axis * (1.0 - kepler.eccentricity) - EARTH_RADIUS;
            if (perigee < 2000) {
                std::cout << "LEO (Low Earth Orbit)";
            } else if (perigee < 35786) {
                std::cout << "MEO (Medium Earth Orbit)";
            } else if (fabs(kepler.inclination) < 1.0) {
                std::cout << "GEO (Geostationary Orbit)";
            } else {
                std::cout << "HEO (High Earth Orbit)";
            }
            std::cout << "\n";
            
        } catch (const std::exception& e) {
            std::cout << "  Error: " << e.what() << "\n";
        }
    }
}

// 交互式程序
void interactiveProgram() {
    std::cout << "TLE to Kepler Elements Converter\n";
    std::cout << "================================\n";
    
    while (true) {
        std::cout << "\nOptions:\n";
        std::cout << "1. Convert TLE from file\n";
        std::cout << "2. Enter TLE manually\n";
        std::cout << "3. Run example conversions\n";
        std::cout << "4. Exit\n";
        std::cout << "Select option: ";
        
        int choice;
        std::cin >> choice;
        std::cin.ignore();  // 清除换行符
        
        switch (choice) {
            case 1: {
                std::cout << "Enter filename: ";
                std::string filename;
                std::getline(std::cin, filename);
                batchConversion(filename);
                break;
            }
            
            case 2: {
                std::string name, line1, line2;
                std::cout << "Enter satellite name (optional): ";
                std::getline(std::cin, name);
                
                std::cout << "Enter TLE line 1: ";
                std::getline(std::cin, line1);
                
                std::cout << "Enter TLE line 2: ";
                std::getline(std::cin, line2);
                
                if (TLE2Kepler::validateTLE(line1, line2)) {
                    try {
                        TLE tle;
                        TLE2Kepler::parseTLE(name, line1, line2, tle);
                        TLE2Kepler::printTLEData(tle);
                        
                        KeplerElements kepler = TLE2Kepler::convertTLEtoKepler(tle);
                        TLE2Kepler::printKeplerElements(kepler);
                    } catch (const std::exception& e) {
                        std::cout << "Error: " << e.what() << "\n";
                    }
                } else {
                    std::cout << "Invalid TLE format\n";
                }
                break;
            }
            
            case 3:
                testTLEConversion();
                break;
                
            case 4:
                std::cout << "Goodbye!\n";
                return;
                
            default:
                std::cout << "Invalid choice. Please try again.\n";
        }
    }
}

int main(int argc, char* argv[]) {
    if (argc > 1) {
        // 命令行模式
        std::string filename = argv[1];
        batchConversion(filename);
    } else {
        // 交互模式
        interactiveProgram();
    }
    
    return 0;
}

9. CMake配置文件

cmake_minimum_required(VERSION 3.10)
project(TLE2Kepler)

set(CMAKE_CXX_STANDARD 11)
set(CMAKE_CXX_STANDARD_REQUIRED ON)

# 设置编译器选项
if(CMAKE_BUILD_TYPE STREQUAL "Release")
    add_compile_options(-O3 -Wall)
else()
    add_compile_options(-O0 -g -Wall)
endif()

# 源文件列表
set(SOURCES
    src/tle2kepler.cpp
    src/kepler_solver.cpp
    src/orbital_calculations.cpp
    src/conversion.cpp
    src/coordinate_conversion.cpp
    src/output_formatter.cpp
    src/julian_date.cpp
    src/main.cpp
)

# 头文件
set(HEADERS
    include/tle2kepler.h
)

# 可执行文件
add_executable(tle2kepler ${SOURCES} ${HEADERS})

# 安装目标
install(TARGETS tle2kepler DESTINATION bin)
install(FILES ${HEADERS} DESTINATION include/tle2kepler)

10. 使用示例和测试数据

10.1 示例TLE文件 (sample.tle)

ISS (ZARYA)
1 25544U 98067A   22001.00000000  .00000000  00000-0  00000-0 0  9999
2 25544  51.6416  33.2520 0006926 359.9521  71.0405 15.49764999300001
GPS BIIR-2  (PRN 13)
1 24876U 97035A   22001.00000000  .00000000  00000-0  00000-0 0  9999
2 24876  55.8093 266.9644 0190932 272.4905  85.8105  2.00562541123456
NOAA 19
1 33591U 09005A   22001.00000000  .00000000  00000-0  00000-0 0  9999
2 33591  99.1762 100.5563 0014088 108.9315 251.3132 14.12404331745634
HST
1 20580U 90037B   22001.00000000  .00000000  00000-0  00000-0 0  9999
2 20580  28.4693 284.8043 0002817 299.7637 123.0579 15.08906432267323

10.2 编译和使用

# 编译
mkdir build
cd build
cmake ..
make

# 运行示例
./tle2kepler

# 处理TLE文件
./tle2kepler ../data/sample.tle

11. 高级功能扩展

11.1 轨道传播器 (orbit_propagator.h)

#ifndef ORBIT_PROPAGATOR_H
#define ORBIT_PROPAGATOR_H

#include "tle2kepler.h"
#include <vector>

class OrbitPropagator {
private:
    KeplerElements initial_elements;
    double initial_epoch;
    
public:
    OrbitPropagator(const KeplerElements& kepler, double epoch)
        : initial_elements(kepler), initial_epoch(epoch) {}
    
    // SGP4/SDP4轨道传播（简化版）
    KeplerElements propagateSGP4(double delta_time_days) {
        KeplerElements propagated = initial_elements;
        
        // 简化传播：只考虑平均运动
        double delta_revs = initial_elements.mean_motion * delta_time_days;
        double new_mean_anomaly = initial_elements.mean_anomaly + 
                                  delta_revs * 360.0;  // 度
        
        // 归一化到[0, 360)
        while (new_mean_anomaly >= 360.0) new_mean_anomaly -= 360.0;
        while (new_mean_anomaly < 0.0) new_mean_anomaly += 360.0;
        
        propagated.mean_anomaly = new_mean_anomaly;
        
        // 重新计算真近点角
        propagated.true_anomaly = TLE2Kepler::meanToTrueAnomaly(
            propagated.mean_anomaly, propagated.eccentricity);
        
        return propagated;
    }
    
    // 生成轨道轨迹
    std::vector<KeplerElements> generateTrajectory(double start_time, 
                                                  double end_time, 
                                                  double time_step) {
        std::vector<KeplerElements> trajectory;
        
        double current_time = start_time;
        while (current_time <= end_time) {
            double delta_days = current_time - initial_epoch;
            KeplerElements current = propagateSGP4(delta_days);
            trajectory.push_back(current);
            
            current_time += time_step;
        }
        
        return trajectory;
    }
};

12. 注意事项

1. 精度限制：此实现是简化版本，实际应用中需要使用SGP4/SDP4等专业轨道模型
2. 坐标系：TLE使用TEME坐标系，而本代码使用简化模型
3. 时间系统：需要正确处理UTC、TAI、TT等时间系统
4. 地球扁率：完整的实现需要考虑J2、J3等地球扁率项
5. 大气阻力：低轨道卫星需要考虑大气阻力影响