Hipparchus-Math / hipparchus

An efficient, general-purpose mathematics components library in the Java programming language
Apache License 2.0
142 stars 41 forks source link

[Feature Request] Kalman Filter - (new module hipparchus-filtering) #32

Closed romgerale closed 6 years ago

romgerale commented 6 years ago

We have been discussing about the necessity of a KalmanFilter in the Hipparchus Luc, I, Maxime).

One alternative, the proposal of this issue, is to create a new submodule called "hipparchus-filtering" that will contain the code from the last version of Apache Commons Math (http://commons.apache.org/proper/commons-math/javadocs/api-3.6/org/apache/commons/math3/filter/package-summary.html) ported to Hipparchus.

A small enhancement will be to change the KalmanFilter implementation from apache because it uses a very strict DEFAULT_RELATIVE_SYMMETRY_THRESHOLD in the "correct" method (namely, in the CholeskyDecomposition). Therefore, it would be great if the KalmanFilter receives, in an additional constructor, the relative symmetry threshold for the CholeskyDecomposition.

The goal of this issue is to allow a centralized discussion between the contributors.

wardev commented 6 years ago

The GaussNewtonOptimizer has a similar issue where it allows the user to choose the matrix decomposition method, but not the tolerances/thresholds. Makes me think there is a use case for a MatrixDecomposer class that encapsulates the algorithm's parameters (e.g. tolerances). Something like:

MatrixDecomposer lu = new LUDecomposer(singularityThreshold);
RealMatrix m = ...;
DecompositionSolver decomposition = lu.decompose(m);
maisonobe commented 6 years ago

We talk with Maxime yesterday at office. I will try to look at the implementation in Orekit next week (starting February 26th) to extract the pure mathematical part from the orbit determination part and integrate it back into Hipparchus. I'll take care to have the various thresholds customizable.

maisonobe commented 6 years ago

I have started work on reimplementing the Kalman. You can find the current Work In Progress in a new branch named kalman. It currently features only an basic implementation of linear kalman, I will add other implementations later on. This corresponds to a new Hipparchus module named filtering.

The API is completely different from the Apache Commons Math one. It mainly consist of two data objects: ProcessEstimate which contains the sate and its covariance, and Measurement which contains the measurement value, covariance and Jacobian. The Kalman estimator will consume a stream of measurements and produce a stream of process estimates (I thought the Java 8 Stream API was really well suited for Kalman).

KalmanEstimator is the base interface and it provides the stream API as a default method, relying on its implementations providing a single step estimate (one measurement produces one new estimate). The first implementation is LinearKalmanEstimator. It is build from a LinearProcess which is an interface users must implement to model their own process (i.e. to provide the matrices used in the linear case, including control). Later on, if this global architecture suits everyone, I will set up an ExtendedKalmanEstimator which will require the user to provide something slightly different (Jacobians matrices for the state, but probably no control at all as I think in non-linear cases the control is in fact embedded in the non-linearity and included in the first Jacobians). This will be this second implementation that we will call from Orekit.

The linear case does compile and even provides results, and there is a set of 3 basic tests (mono-dimensional, no control, estimation of a random constant). I don't know yet if the results are correct. I'm not really happy with the global behavior, so I have to look at it more closely, add news tests and compare with other implementations.

I would be happy if you could take a look ot this and comment.

maisonobe commented 6 years ago

Another comment: I have added the MatrixDecomposer feature as proposed by Evan, and uses it here, so matrix decomposition can indeed be configured by user.

maisonobe commented 6 years ago

I have validated the linear Kalman filter against the Apache Commons Math version. Differences between the two implementations at at double precision noise level, despite the they have absolutely nothing in common. Now I will start implementing a non-linear version.

maisonobe commented 6 years ago

Here is what I have in mind to handle both linear and non-linear Kalman. kalman-design

maisonobe commented 6 years ago

The kalman branch has been merged into master. Everything seems to work properly. The final class diagram is only slightly different from the previous one, mainly in the non-linear case. The change allowed to let the user postpone the decision to reject the measurement to between the prediction and correction steps, when the innovation covariance matrix is available. A typical use case is that we reject outlier measurements in orbit determination based on the ratio between residuals and covariance, so we need to have a correct estimation of the innovation covariance matrix which is only available after the prediction step. kalman