Hello, I would like to use orekit to compare some propagators. I have data from telescope and TLE frames before and after those. My scheme would be:
- change TLE to position, velocity vector with TLEPropagator
- use data with kalman filter and different propagators
- propagate it to data of next avaible TLE frame
- compare it with second changed into vector TLE frame.
I wrote this
import orekit
import math
import matplotlib.pyplot as plt
import numpy as np
from orekit.pyhelpers import setup_orekit_curdir
from org.orekit.frames import FramesFactory
from org.orekit.time import TimeScalesFactory, AbsoluteDate, DateComponents, TimeComponents
from org.orekit.propagation.analytical.tle import TLE, TLEPropagator
from org.orekit.propagation.conversion import NumericalPropagatorBuilder, DormandPrince853IntegratorBuilder
from org.orekit.orbits import PositionAngleType
from org.hipparchus.linear import Array2DRowRealMatrix
from org.orekit.estimation.sequential import KalmanEstimatorBuilder, ConstantProcessNoise
from org.orekit.estimation.measurements import AngularRaDec, ObservableSatellite, GroundStation
from org.orekit.bodies import OneAxisEllipsoid, GeodeticPoint
from org.orekit.frames import TopocentricFrame
from org.orekit.utils import IERSConventions, Constants
from org.hipparchus.geometry.euclidean.threed import Vector3D
from org.orekit.utils import PVCoordinates
from org.orekit.propagation.analytical import KeplerianPropagator
from org.orekit.orbits import CartesianOrbit
from math import atan2
orekit.initVM()
setup_orekit_curdir()
line1 = "1 44542U 00000A 24141.59902293 .00000000 00000-0 00000-0 0 9992"
line2 = "2 44542 54.4447 208.2914 0007327 347.2490 12.8128 1.86231282 00"
tle = TLE(line1, line2)
tle_propagator = TLEPropagator.selectExtrapolator(tle)
frame = FramesFactory.getEME2000()
initial_state = tle_propagator.getInitialState()
initial_position = initial_state.getPVCoordinates().getPosition()
initial_velocity = initial_state.getPVCoordinates().getVelocity()
print("Wektor stanu po przeliczeniu z TLE:")
print(f"Pozycja (m): X={initial_position.getX():.3f}, Y={initial_position.getY():.3f}, Z={initial_position.getZ():.3f}")
print(f"Prędkość (m/s): VX={initial_velocity.getX():.3f}, VY={initial_velocity.getY():.3f}, VZ={initial_velocity.getZ():.3f}")
earth = OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
Constants.WGS84_EARTH_FLATTENING,
FramesFactory.getITRF(IERSConventions.IERS_2010, True))
longitude = math.radians(149.069793)
latitude = math.radians(-31.27344)
altitude = 1115.0
station_position = GeodeticPoint(latitude, longitude, altitude)
station_frame = TopocentricFrame(earth, station_position, "Esrange")
ground_station = GroundStation(station_frame)
min_step = 0.001
max_step = 300.0
position_tolerance = 10.0
integrator_builder = DormandPrince853IntegratorBuilder(min_step, max_step, position_tolerance)
propagator_builder = NumericalPropagatorBuilder(
initial_state.getOrbit(), integrator_builder, PositionAngleType.MEAN, position_tolerance
)
covariance_matrix = Array2DRowRealMatrix(6, 6)
for i in range(6):
covariance_matrix.setEntry(i, i, 1e-5)
process_noise = ConstantProcessNoise(covariance_matrix, covariance_matrix)
kalman_builder = KalmanEstimatorBuilder()
kalman_builder.addPropagationConfiguration(propagator_builder, process_noise)
kalman_estimator = kalman_builder.build()
observations = [
("2024-05-20T14:45:44.752829", 206.44878323599582, 17.776891936217453),
("2024-05-20T14:45:45.553957", 206.45218390619834, 17.782686863708093),
("2024-05-20T14:45:46.355257", 206.45575566125694, 17.788097557353804),
("2024-05-20T14:45:47.156427", 206.45959448811857, 17.793920431927166),
("2024-05-20T14:45:47.957615", 206.46334901411694, 17.799551151912503),
("2024-05-20T14:45:48.758819", 206.4663549931162, 17.805271879450217),
("2024-05-20T14:45:49.560083", 206.47046313292233, 17.810888131897244),
("2024-05-20T14:45:50.361241", 206.47370698057398, 17.816519271633158),
("2024-05-20T14:45:51.162418", 206.47810090124472, 17.822224086108644),
("2024-05-20T14:45:51.963735", 206.4812888168711, 17.82785727176931),
("2024-05-20T14:45:52.764924", 206.48504044130928, 17.83357182465419),
("2024-05-20T14:45:53.566107", 206.4882699636153, 17.839390706585426),
]
def parse_date(date_str, utc):
date_part, time_part = date_str.split("T")
year, month, day = map(int, date_part.split("-"))
hour, minute, second = map(float, time_part.split(":"))
total_seconds = hour * 3600 + minute * 60 + second # Sekundy dnia
return AbsoluteDate(DateComponents(year, month, day), TimeComponents(total_seconds), utc)
observable_satellite = ObservableSatellite(0)
updated_states = []
for observation in observations:
utc = TimeScalesFactory.getUTC()
date_str, ra, dec = observation
date = parse_date(date_str, utc)
ra_rad = math.radians(ra)
dec_rad = math.radians(dec)
angular_measurement = AngularRaDec(
ground_station, frame, date, [ra_rad, dec_rad], [1.0, 1.0], [1.0, 1.0], observable_satellite
)
kalman_estimator.estimationStep(angular_measurement)
estimated_state = kalman_estimator.getPhysicalEstimatedState()
position = estimated_state.getSubVector(0, 3) # Pozycja: X, Y, Z
velocity = estimated_state.getSubVector(3, 3) # Prędkość: VX, VY, VZ
updated_states.append([position, velocity])
print(f"\nZaktualizowany stan dla obserwacji z {date}:")
print(f"Pozycja (m): X={position.getEntry(0):.3f}, Y={position.getEntry(1):.3f}, Z={position.getEntry(2):.3f}")
print(f"Prędkość (m/s): VX={velocity.getEntry(0):.3f}, VY={velocity.getEntry(1):.3f}, VZ={velocity.getEntry(2):.3f}")
target_date = AbsoluteDate(
DateComponents(2024, 5, 21),
TimeComponents(3 * 3600 + 16 * 60 + 29.0),
TimeScalesFactory.getUTC()
)
propagator = propagator_builder.buildPropagator(propagator_builder.getSelectedNormalizedParameters())
propagated_state = propagator.propagate(target_date)
propagated_position = propagated_state.getPVCoordinates().getPosition()
propagated_velocity = propagated_state.getPVCoordinates().getVelocity()
print("\nStan satelity na dzień 2024-05-21 03:16:29 UTC:")
print(f"Pozycja (m): X={propagated_position.getX():.3f}, Y={propagated_position.getY():.3f}, Z={propagated_position.getZ():.3f}")
print(f"Prędkość (m/s): VX={propagated_velocity.getX():.3f}, VY={propagated_velocity.getY():.3f}, VZ={propagated_velocity.getZ():.3f}")
`
NumericalPropagator works (I guess) but I can’t compile other, I tried Keplerian and DSST without succes
Is it possible? I would like to use other also 
---------------------------------------------------------------------------
InvalidArgsError Traceback (most recent call last)
Cell In[1], line 110
108 process_noise = ConstantProcessNoise(covariance_matrix, covariance_matrix)
109 kalman_builder = KalmanEstimatorBuilder()
--> 110 kalman_builder.addPropagationConfiguration(keplerian_propagator, process_noise)
111 kalman_estimator = kalman_builder.build()
113 # Obserwacje
InvalidArgsError: (<class 'org.orekit.estimation.sequential.KalmanEstimatorBuilder'>, 'addPropagationConfiguration', (<KeplerianPropagator: org.orekit.propagation.analytical.KeplerianPropagator@3701e6e4>, <ConstantProcessNoise: org.orekit.estimation.sequential.ConstantProcessNoise@37ce3644>))
Could You evaluate my plan and help me with propagators?