Setting Kalman filter

Orekit Python Wrapper usage
1

Hello, I’m trying to estimate orbit using KalmanEstimator, but when trying to add measurement i get error like this

---------------------------------------------------------------------------
JavaError                                 Traceback (most recent call last)
Cell In[1], line 54
     47 measurements = [
     48     Position(initial_date.shiftedBy(300.0), Vector3D(7100000.0, 10000.0, 20000.0), position_error, 1.0, None),
     49     Position(initial_date.shiftedBy(600.0), Vector3D(7200000.0, 20000.0, 30000.0), position_error, 1.0, None),
     50     Position(initial_date.shiftedBy(900.0), Vector3D(7300000.0, 30000.0, 40000.0), position_error, 1.0, None),
     51 ]
     53 for measurement in measurements:
---> 54     kalman_estimator.estimationStep(measurement)
     56 estimated_orbit = kalman_estimator.getEstimatedPropagator().getInitialState().getOrbit()
     58 print(estimated_orbit.getA())

JavaError: <super: <class 'JavaError'>, <JavaError object>>
    Java stacktrace:
java.lang.NullPointerException
	at org.orekit.estimation.sequential.KalmanEstimatorUtil.filterRelevant(KalmanEstimatorUtil.java:199)
	at org.orekit.estimation.sequential.KalmanModel.getEvolution(KalmanModel.java:713)
	at org.orekit.estimation.sequential.KalmanModel.getEvolution(KalmanModel.java:51)
	at org.hipparchus.filtering.kalman.extended.ExtendedKalmanFilter.estimationStep(ExtendedKalmanFilter.java:56)
	at org.orekit.estimation.sequential.KalmanEstimator.estimationStep(KalmanEstimator.java:142)

Could you help me with it?

import orekit
from orekit.pyhelpers import setup_orekit_curdir
from org.orekit.propagation.conversion import DormandPrince853IntegratorBuilder, NumericalPropagatorBuilder
from org.orekit.estimation.sequential import KalmanEstimatorBuilder, ConstantProcessNoise
from org.orekit.orbits import KeplerianOrbit, PositionAngleType
from org.orekit.utils import Constants, PVCoordinates
from org.orekit.time import AbsoluteDate, TimeScalesFactory
from org.orekit.frames import FramesFactory
from org.hipparchus.linear import QRDecomposer, BlockRealMatrix
from org.orekit.estimation.measurements import Position
from org.hipparchus.geometry.euclidean.threed import Vector3D

orekit.initVM()
setup_orekit_curdir()

mu = Constants.WGS84_EARTH_MU
frame = FramesFactory.getEME2000()

# Kontekst czasu
utc = TimeScalesFactory.getUTC()
initial_date = AbsoluteDate(2024, 4, 25, 20, 52, 0.0, utc)

position = Vector3D(-5915780.397691450, -966818.290504886, -3206121.242323730)
velocity = Vector3D(-1651.90009192014, -5738.13475197905, 4790.80932769098)

pv_coordinates = PVCoordinates(position, velocity)

k_orbit = KeplerianOrbit(pv_coordinates, frame, initial_date, mu)

min_step = 0.01
max_step = 600.0
position_error = 0.001

integrator_builder = DormandPrince853IntegratorBuilder(min_step, max_step, position_error)

propagator_builder = NumericalPropagatorBuilder(k_orbit, integrator_builder, PositionAngleType.TRUE, 1.0)

covariance_matrix = BlockRealMatrix(6, 6)
constant_process_noise = ConstantProcessNoise(covariance_matrix, covariance_matrix)

kalman_builder = KalmanEstimatorBuilder()
kalman_builder.addPropagationConfiguration(propagator_builder, constant_process_noise)
kalman_builder.decomposer(QRDecomposer(1.0e-15))

kalman_estimator = kalman_builder.build()

measurements = [
    Position(initial_date.shiftedBy(300.0), Vector3D(7100000.0, 10000.0, 20000.0), position_error, 1.0, None),
    Position(initial_date.shiftedBy(600.0), Vector3D(7200000.0, 20000.0, 30000.0), position_error, 1.0, None),
    Position(initial_date.shiftedBy(900.0), Vector3D(7300000.0, 30000.0, 40000.0), position_error, 1.0, None),
]

for measurement in measurements:
    kalman_estimator.estimationStep(measurement)

estimated_orbit = kalman_estimator.getEstimatedPropagator().getInitialState().getOrbit()

print(estimated_orbit.getA())
2

Hi @Tomuch,

Welcome to Orekit’s forum!

I think your problem is that the ObservableSatellite is None in your measurements.

Replace:

measurements = [
     Position(initial_date.shiftedBy(300.0), Vector3D(7100000.0, 10000.0, 20000.0), position_error, 1.0, None),
     Position(initial_date.shiftedBy(600.0), Vector3D(7200000.0, 20000.0, 30000.0), position_error, 1.0, None),
     Position(initial_date.shiftedBy(900.0), Vector3D(7300000.0, 30000.0, 40000.0), position_error, 1.0, None),
 ]

With:

from org.orekit.estimation.measurements import ObservableSatellite
my_sat = ObservableSatellite(0)

measurements = [
     Position(initial_date.shiftedBy(300.0), Vector3D(7100000.0, 10000.0, 20000.0), position_error, 1.0, my_sat),
     Position(initial_date.shiftedBy(600.0), Vector3D(7200000.0, 20000.0, 30000.0), position_error, 1.0, my_sat),
     Position(initial_date.shiftedBy(900.0), Vector3D(7300000.0, 30000.0, 40000.0), position_error, 1.0, my_sat),
 ]

Hope this helps,
Maxime

3

That solved my problem, thanks

Could you tell me what it is actually? Is it just some kind of name or ID?

4

Yes it’s a number (or ID) used to assign a measurement to a satellite in case of multi-satellites orbit determination.
Maybe we should add the possibility to assign a name to it by the way

1 Like