#!/usr/bin/env python3
"""
                                tle2ssc
    Convert TLE orbit to Solar System Catalog (SSC) format for Celestia.

    Copyright (C) 2022, 2023, Jeff Moe
    Authors: Jeff Moe

    This program is free software: you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 2 of the License, or
    (at your option) any later version.

    This program is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with this program. If not, see <https://www.gnu.org/licenses/>.
"""

# Usage:
# tle2ssc
# Example:
# ./tle2ssc > /usr/share/celestia/extras-standard/gnss/gnss.ssc

import os
import math
from math import remainder
from datetime import datetime
from sgp4.api import Satrec, WGS72
from skyfield.api import EarthSatellite, load, wgs84

xpdotp = 1440.0 / (2.0 * math.pi)

# Use GNSS to get all four systems
satellites_url = 'https://celestrak.com/NORAD/elements/gnss.txt'
#satellites_url = 'https://celestrak.com/NORAD/elements/galileo.txt'
#satellites_url = 'https://celestrak.com/NORAD/elements/gps-ops.txt'
#satellites_url = 'https://celestrak.com/NORAD/elements/glo-ops.txt'
#satellites_url = 'https://celestrak.com/NORAD/elements/beidou.txt'

satellites = load.tle_file(satellites_url,reload=True)

ts = load.timescale()
t = ts.now()

for satellite in satellites:
	satellite_name=satellite.name

# Static radius
	satellite_radius=0.005

# Epoch XXX
	#satellite_epoch=satellite.model.jdsatepoch
	satellite_epoch=t.tdb

#    The unique satellite NORAD catalog number given in the TLE file.
	satellite_number=satellite.model.satnum

#    Satellite classification, or else 'U' for “Unknown”
	satellite_classification=satellite.model.classification

#    International designator
	satellite_intldesg=satellite.model.intldesg

#    Full four-digit year of this element set’s epoch moment.
	satellite_epochyr=satellite.model.epochyr

#    Fractional days into the year of the epoch moment.
	satellite_epochdays=satellite.model.epochdays

#    Julian date of the epoch (computed from epochyr and epochdays).
	satellite_jdsatepoch=satellite.model.jdsatepoch

#    First time derivative of the mean motion (ignored by SGP4).
	satellite_ndot=satellite.model.ndot

#    Second time derivative of the mean motion (ignored by SGP4).
	satellite_nddot=satellite.model.nddot

#    Ballistic drag coefficient B* in inverse earth radii.
	satellite_bstar=satellite.model.bstar

#    Ephemeris type (ignored by SGP4 as determination now automatic)
	satellite_ephtype=satellite.model.ephtype

#    Element number
	satellite_elnum=satellite.model.elnum

#    Inclination in radians. Convert radians to degrees.
	satellite_inclination=math.degrees(satellite.model.inclo)

# Satellite Inclination and Obliquity are the same
	satellite_obliquity=math.degrees(satellite.model.inclo)

#    Right ascension of ascending node in radians. Convert to degrees.
	satellite_ascending_node=math.degrees(satellite.model.nodeo)

# Equator Ascending Node, in degrees
	satellite_equator_ascending_node=math.degrees(satellite.model.nodeo)

#    Eccentricity.
	satellite_eccentricity=satellite.model.ecco

#    Argument of perigee in radians.
	satellite_arg_of_pericenter=math.degrees(satellite.model.argpo)

#    Mean anomaly in radians.
	satellite_mean_anomaly=math.degrees(satellite.model.mo)

#    Mean motion in radians per minute.
	satellite_period=(1 / (satellite.model.no_kozai * xpdotp))

#    Revolution number at epoch [Revs]
	satellite_revnum=satellite.model.revnum

# Semi-Major Axis
	satellite_semimajor_axis=pow((pow(satellite_period,2) * 75371000000000),0.33333333)

# Rotation Offset
	satellite_rotation_offset=((satellite_arg_of_pericenter)+(satellite_mean_anomaly)+360*
		((2451545-satellite_epochdays)/satellite_period-
		round((2451545-satellite_epochdays)/satellite_period)))

# Output for .ssc file
	print('"', satellite_name, '" "Sol/Earth" {',sep="")
	print('    Class "spacecraft"')
	print('    Mesh "galileo-gnss.cmod"')
	print('    Radius', "%.3f" %satellite_radius)
	print()
	print('  EllipticalOrbit {')
	print('    Epoch', "%.8f" %satellite_epoch)
	print('    Period', "%.8f" %satellite_period)
	print('    SemiMajorAxis', "%.3f" %satellite_semimajor_axis)
	print('    Eccentricity', "%.8f" %satellite_eccentricity)
	print('    Inclination', "%.4f" %satellite_inclination)
	print('    AscendingNode', "%.4f" %satellite_ascending_node)
	print('    ArgOfPericenter', "%.4f" %satellite_arg_of_pericenter)
	print('    MeanAnomaly', "%.4f" %satellite_mean_anomaly)
	print('  }')
	print('  Obliquity', "%.4f" %satellite_obliquity)
	print('  EquatorAscendingNode', "%.4f" %satellite_equator_ascending_node)
	print('  RotationOffset', "%.4f" %satellite_rotation_offset)
	print('  # Orientation  [  ]')
	print('}')
	print()