Skip to content

GitLab

Commit cd833792 authored by Bryan Cazabonne's avatar Bryan Cazabonne
Browse files

Added unit tests for orbit determination based on Brouwer-Lyddane.

parent 1470b4fd
Expand all Hide whitespace changes
Inline Side-by-side
src/test/java/org/orekit/estimation/BrouwerLyddaneContext.java 0 → 100644
View file @ cd833792
package org.orekit.estimation;
import java.util.List;
import java.util.Map;
import org.hipparchus.geometry.euclidean.threed.Vector3D;
import org.hipparchus.util.FastMath;
import org.orekit.bodies.CelestialBody;
import org.orekit.bodies.GeodeticPoint;
import org.orekit.bodies.OneAxisEllipsoid;
import org.orekit.estimation.measurements.GroundStation;
import org.orekit.forces.gravity.potential.UnnormalizedSphericalHarmonicsProvider;
import org.orekit.frames.TopocentricFrame;
import org.orekit.models.earth.displacement.StationDisplacement;
import org.orekit.orbits.CartesianOrbit;
import org.orekit.orbits.Orbit;
import org.orekit.orbits.PositionAngle;
import org.orekit.propagation.analytical.BrouwerLyddanePropagator;
import org.orekit.propagation.conversion.BrouwerLyddanePropagatorBuilder;
import org.orekit.time.TimeScale;
import org.orekit.time.UT1Scale;
import org.orekit.utils.IERSConventions;
import org.orekit.utils.PVCoordinates;
public class BrouwerLyddaneContext implements StationDataProvider {
public IERSConventions conventions;
public OneAxisEllipsoid earth;
public CelestialBody sun;
public CelestialBody moon;
public UnnormalizedSphericalHarmonicsProvider gravity;
public TimeScale utc;
public UT1Scale ut1;
public Orbit initialOrbit;
public StationDisplacement[] displacements;
public List<GroundStation> stations;
// Stations for turn-around range
// Map entry = primary station
// Map value = secondary station associated
public Map<GroundStation, GroundStation> TARstations;
public BrouwerLyddanePropagatorBuilder createBuilder(final PositionAngle angleType, final boolean perfectStart, final double dP) {
final Orbit startOrbit;
if (perfectStart) {
// orbit estimation will start from a perfect orbit
startOrbit = initialOrbit;
} else {
// orbit estimation will start from a wrong point
final Vector3D initialPosition = initialOrbit.getPVCoordinates().getPosition();
final Vector3D initialVelocity = initialOrbit.getPVCoordinates().getVelocity();
final Vector3D wrongPosition = initialPosition.add(new Vector3D(1000.0, 0, 0));
final Vector3D wrongVelocity = initialVelocity.add(new Vector3D(0, 0, 0.01));
startOrbit = new CartesianOrbit(new PVCoordinates(wrongPosition, wrongVelocity),
initialOrbit.getFrame(),
initialOrbit.getDate(),
initialOrbit.getMu());
}
// Initialize builder
final BrouwerLyddanePropagatorBuilder propagatorBuilder =
new BrouwerLyddanePropagatorBuilder(startOrbit, gravity, angleType, dP, BrouwerLyddanePropagator.M2);
// Return
return propagatorBuilder;
}
GroundStation createStation(double latitudeInDegrees, double longitudeInDegrees,
double altitude, String name) {
final GeodeticPoint gp = new GeodeticPoint(FastMath.toRadians(latitudeInDegrees),
FastMath.toRadians(longitudeInDegrees),
altitude);
return new GroundStation(new TopocentricFrame(earth, gp, name),
ut1.getEOPHistory(), displacements);
}
@Override
public List<GroundStation> getStations() {
return stations;
}
}
src/test/java/org/orekit/estimation/leastsquares/BrouwerLyddaneBatchLSEstimatorTest.java 0 → 100644
View file @ cd833792
/* Copyright 2002-2022 CS GROUP
* Licensed to CS GROUP (CS) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* CS licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.orekit.estimation.leastsquares;
import java.util.ArrayList;
import java.util.List;
import org.hipparchus.geometry.euclidean.threed.Vector3D;
import org.hipparchus.linear.RealMatrix;
import org.hipparchus.optim.nonlinear.vector.leastsquares.LeastSquaresProblem.Evaluation;
import org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer;
import org.junit.Assert;
import org.junit.Test;
import org.orekit.attitudes.LofOffset;
import org.orekit.errors.OrekitException;
import org.orekit.errors.OrekitMessages;
import org.orekit.estimation.BrouwerLyddaneContext;
import org.orekit.estimation.BrouwerLyddaneEstimationTestUtils;
import org.orekit.estimation.measurements.EstimationsProvider;
import org.orekit.estimation.measurements.GroundStation;
import org.orekit.estimation.measurements.ObservedMeasurement;
import org.orekit.estimation.measurements.PVMeasurementCreator;
import org.orekit.estimation.measurements.Range;
import org.orekit.estimation.measurements.RangeMeasurementCreator;
import org.orekit.estimation.measurements.RangeRateMeasurementCreator;
import org.orekit.estimation.measurements.modifiers.OnBoardAntennaRangeModifier;
import org.orekit.frames.LOFType;
import org.orekit.orbits.Orbit;
import org.orekit.orbits.PositionAngle;
import org.orekit.propagation.Propagator;
import org.orekit.propagation.conversion.BrouwerLyddanePropagatorBuilder;
import org.orekit.utils.ParameterDriversList;
public class BrouwerLyddaneBatchLSEstimatorTest {
/**
* Perfect PV measurements with a perfect start
*/
@Test
public void testPV() {
BrouwerLyddaneContext context = BrouwerLyddaneEstimationTestUtils.eccentricContext("regular-data:potential:tides");
final BrouwerLyddanePropagatorBuilder propagatorBuilder = context.createBuilder(PositionAngle.MEAN, true, 1.0);
// create perfect PV measurements
final Propagator propagator = BrouwerLyddaneEstimationTestUtils.createPropagator(context.initialOrbit,
propagatorBuilder);
final List<ObservedMeasurement<?>> measurements =
BrouwerLyddaneEstimationTestUtils.createMeasurements(propagator,
new PVMeasurementCreator(),
0.0, 1.0, 300.0);
// create orbit estimator
final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
propagatorBuilder);
for (final ObservedMeasurement<?> measurement : measurements) {
estimator.addMeasurement(measurement);
}
estimator.setParametersConvergenceThreshold(1.0e-2);
estimator.setMaxIterations(10);
estimator.setMaxEvaluations(20);
BrouwerLyddaneEstimationTestUtils.checkFit(context, estimator, 1, 1,
0.0, 1.0e-15,
0.0, 1.0e-15,
0.0, 1.0e-15,
0.0, 1.0e-15);
RealMatrix normalizedCovariances = estimator.getOptimum().getCovariances(1.0e-10);
RealMatrix physicalCovariances = estimator.getPhysicalCovariances(1.0e-10);
Assert.assertEquals(6, normalizedCovariances.getRowDimension());
Assert.assertEquals(6, normalizedCovariances.getColumnDimension());
Assert.assertEquals(6, physicalCovariances.getRowDimension());
Assert.assertEquals(6, physicalCovariances.getColumnDimension());
Assert.assertEquals(0.0, physicalCovariances.getEntry(0, 0), 1.7e-15);
}
/** Test PV measurements generation and backward propagation in least-square orbit determination. */
@Test
public void testKeplerPVBackward() {
BrouwerLyddaneContext context = BrouwerLyddaneEstimationTestUtils.eccentricContext("regular-data:potential:tides");
final BrouwerLyddanePropagatorBuilder propagatorBuilder =
context.createBuilder(PositionAngle.MEAN, true, 1.0);
// create perfect PV measurements
final Propagator propagator = BrouwerLyddaneEstimationTestUtils.createPropagator(context.initialOrbit,
propagatorBuilder);
final List<ObservedMeasurement<?>> measurements =
BrouwerLyddaneEstimationTestUtils.createMeasurements(propagator,
new PVMeasurementCreator(),
0.0, -1.0, 300.0);
// create orbit estimator
final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
propagatorBuilder);
for (final ObservedMeasurement<?> measurement : measurements) {
estimator.addMeasurement(measurement);
}
estimator.setParametersConvergenceThreshold(1.0e-2);
estimator.setMaxIterations(10);
estimator.setMaxEvaluations(20);
BrouwerLyddaneEstimationTestUtils.checkFit(context, estimator, 1, 1,
0.0, 1.0e-15,
0.0, 1.0e-15,
0.0, 1.0e-15,
0.0, 1.0e-15);
RealMatrix normalizedCovariances = estimator.getOptimum().getCovariances(1.0e-10);
RealMatrix physicalCovariances = estimator.getPhysicalCovariances(1.0e-10);
Assert.assertEquals(6, normalizedCovariances.getRowDimension());
Assert.assertEquals(6, normalizedCovariances.getColumnDimension());
Assert.assertEquals(6, physicalCovariances.getRowDimension());
Assert.assertEquals(6, physicalCovariances.getColumnDimension());
Assert.assertEquals(0.0, physicalCovariances.getEntry(0, 0), 1.7e-15);
}
/**
* Perfect range measurements with a perfect start
*/
@Test
public void testKeplerRange() {
BrouwerLyddaneContext context = BrouwerLyddaneEstimationTestUtils.eccentricContext("regular-data:potential:tides");
final BrouwerLyddanePropagatorBuilder propagatorBuilder =
context.createBuilder(PositionAngle.MEAN, true, 1.0);
// create perfect range measurements
final Propagator propagator = BrouwerLyddaneEstimationTestUtils.createPropagator(context.initialOrbit,
propagatorBuilder);
final List<ObservedMeasurement<?>> measurements =
BrouwerLyddaneEstimationTestUtils.createMeasurements(propagator,
new RangeMeasurementCreator(context),
1.0, 3.0, 300.0);
// create orbit estimator
final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
propagatorBuilder);
for (final ObservedMeasurement<?> range : measurements) {
estimator.addMeasurement(range);
}
estimator.setParametersConvergenceThreshold(1.0e-2);
estimator.setMaxIterations(10);
estimator.setMaxEvaluations(20);
BrouwerLyddaneEstimationTestUtils.checkFit(context, estimator, 1, 2,
0.0, 1.1e-4,
0.0, 1.8e-4,
0.0, 1.4e-5,
0.0, 1.3e-8);
}
/**
* Perfect range measurements with a perfect start and an on-board antenna range offset
*/
@Test
public void testKeplerRangeWithOnBoardAntennaOffset() {
BrouwerLyddaneContext context = BrouwerLyddaneEstimationTestUtils.eccentricContext("regular-data:potential:tides");
final BrouwerLyddanePropagatorBuilder propagatorBuilder =
context.createBuilder(PositionAngle.MEAN, true, 1.0);
propagatorBuilder.setAttitudeProvider(new LofOffset(propagatorBuilder.getFrame(), LOFType.LVLH));
final Vector3D antennaPhaseCenter = new Vector3D(-1.2, 2.3, -0.7);
// create perfect range measurements with antenna offset
final Propagator propagator = BrouwerLyddaneEstimationTestUtils.createPropagator(context.initialOrbit,
propagatorBuilder);
final List<ObservedMeasurement<?>> measurements =
BrouwerLyddaneEstimationTestUtils.createMeasurements(propagator,
new RangeMeasurementCreator(context, antennaPhaseCenter),
1.0, 3.0, 300.0);
// create orbit estimator
final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
propagatorBuilder);
final OnBoardAntennaRangeModifier obaModifier = new OnBoardAntennaRangeModifier(antennaPhaseCenter);
for (final ObservedMeasurement<?> range : measurements) {
((Range) range).addModifier(obaModifier);
estimator.addMeasurement(range);
}
estimator.setParametersConvergenceThreshold(1.0e-2);
estimator.setMaxIterations(10);
estimator.setMaxEvaluations(20);
BrouwerLyddaneEstimationTestUtils.checkFit(context, estimator, 1, 2,
0.0, 1.2e-4,
0.0, 2.6e-4,
0.0, 1.4e-5,
0.0, 1.3e-8);
}
/**
* Perfect range rate measurements with a perfect start
*/
@Test
public void testKeplerRangeRate() {
BrouwerLyddaneContext context = BrouwerLyddaneEstimationTestUtils.eccentricContext("regular-data:potential:tides");
final BrouwerLyddanePropagatorBuilder propagatorBuilder =
context.createBuilder(PositionAngle.MEAN, true, 1.0);
// create perfect range rate measurements
final Propagator propagator = BrouwerLyddaneEstimationTestUtils.createPropagator(context.initialOrbit,
propagatorBuilder);
final double groundClockDrift = 4.8e-9;
for (final GroundStation station : context.stations) {
station.getClockDriftDriver().setValue(groundClockDrift);
}
final double satClkDrift = 3.2e-10;
final List<ObservedMeasurement<?>> measurements1 =
BrouwerLyddaneEstimationTestUtils.createMeasurements(propagator,
new RangeRateMeasurementCreator(context, false, satClkDrift),
1.0, 3.0, 300.0);
final List<ObservedMeasurement<?>> measurements = new ArrayList<ObservedMeasurement<?>>();
measurements.addAll(measurements1);
// create orbit estimator
final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
propagatorBuilder);
for (final ObservedMeasurement<?> rangerate : measurements) {
estimator.addMeasurement(rangerate);
}
estimator.setParametersConvergenceThreshold(1.0e-3);
estimator.setMaxIterations(10);
estimator.setMaxEvaluations(20);
BrouwerLyddaneEstimationTestUtils.checkFit(context, estimator, 1, 2,
0.0, 7.9e-8,
0.0, 1.1e-7,
0.0, 1.4e-5,
0.0, 1.4e-8);
}
@Test
public void testWrappedException() {
BrouwerLyddaneContext context = BrouwerLyddaneEstimationTestUtils.eccentricContext("regular-data:potential:tides");
final BrouwerLyddanePropagatorBuilder propagatorBuilder =
context.createBuilder(PositionAngle.MEAN, true, 1.0);
// create perfect range measurements
final Propagator propagator = BrouwerLyddaneEstimationTestUtils.createPropagator(context.initialOrbit,
propagatorBuilder);
final List<ObservedMeasurement<?>> measurements =
BrouwerLyddaneEstimationTestUtils.createMeasurements(propagator,
new RangeMeasurementCreator(context),
1.0, 3.0, 300.0);
// create orbit estimator
final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
propagatorBuilder);
for (final ObservedMeasurement<?> range : measurements) {
estimator.addMeasurement(range);
}
estimator.setParametersConvergenceThreshold(1.0e-2);
estimator.setMaxIterations(10);
estimator.setMaxEvaluations(20);
estimator.setObserver(new BatchLSObserver() {
/** {@inheritDoc} */
@Override
public void evaluationPerformed(int iterationsCount, int evaluationscount,
Orbit[] orbits,
ParameterDriversList estimatedOrbitalParameters,
ParameterDriversList estimatedPropagatorParameters,
ParameterDriversList estimatedMeasurementsParameters,
EstimationsProvider evaluationsProvider, Evaluation lspEvaluation) {
throw new DummyException();
}
});
try {
BrouwerLyddaneEstimationTestUtils.checkFit(context, estimator, 3, 4,
0.0, 1.5e-6,
0.0, 3.2e-6,
0.0, 3.8e-7,
0.0, 1.5e-10);
Assert.fail("an exception should have been thrown");
} catch (DummyException de) {
// expected
}
}
private static class DummyException extends OrekitException {
private static final long serialVersionUID = 1L;
public DummyException() {
super(OrekitMessages.INTERNAL_ERROR);
}
}
}
src/test/java/org/orekit/estimation/sequential/BrouwerLyddaneKalmanEstimatorTest.java 0 → 100644
View file @ cd833792
package org.orekit.estimation.sequential;
import java.util.List;
import org.hipparchus.linear.MatrixUtils;
import org.hipparchus.linear.RealMatrix;
import org.junit.Test;
import org.orekit.estimation.BrouwerLyddaneContext;
import org.orekit.estimation.BrouwerLyddaneEstimationTestUtils;
import org.orekit.estimation.measurements.ObservedMeasurement;
import org.orekit.estimation.measurements.PVMeasurementCreator;
import org.orekit.orbits.Orbit;
import org.orekit.orbits.PositionAngle;
import org.orekit.propagation.Propagator;
import org.orekit.propagation.analytical.BrouwerLyddanePropagator;
import org.orekit.propagation.conversion.BrouwerLyddanePropagatorBuilder;
public class BrouwerLyddaneKalmanEstimatorTest {
/**
* Perfect PV measurements with a perfect start
* Keplerian formalism
*/
@Test
public void testKeplerianPV() {
// Create context
BrouwerLyddaneContext context = BrouwerLyddaneEstimationTestUtils.eccentricContext("regular-data:potential:tides");
// Create initial orbit and propagator builder
final PositionAngle positionAngle = PositionAngle.TRUE;
final boolean perfectStart = true;
final double dP = 1.;
final BrouwerLyddanePropagatorBuilder propagatorBuilder =
context.createBuilder(positionAngle, perfectStart, dP);
// Create perfect PV measurements
final Propagator propagator = BrouwerLyddaneEstimationTestUtils.createPropagator(context.initialOrbit,
propagatorBuilder);
final List<ObservedMeasurement<?>> measurements =
BrouwerLyddaneEstimationTestUtils.createMeasurements(propagator,
new PVMeasurementCreator(),
0.0, 3.0, 300.0);
// Reference propagator for estimation performances
final BrouwerLyddanePropagator referencePropagator = propagatorBuilder.
buildPropagator(propagatorBuilder.getSelectedNormalizedParameters());
// Reference position/velocity at last measurement date
final Orbit refOrbit = referencePropagator.
propagate(measurements.get(measurements.size()-1).getDate()).getOrbit();
// Covariance matrix initialization
final RealMatrix initialP = MatrixUtils.createRealDiagonalMatrix(new double [] {
1e-2, 1e-2, 1e-2, 1e-5, 1e-5, 1e-5
});
// Process noise matrix
RealMatrix Q = MatrixUtils.createRealDiagonalMatrix(new double [] {
1.e-8, 1.e-8, 1.e-8, 1.e-8, 1.e-8, 1.e-8
});
// Build the Kalman filter
final KalmanEstimator kalman = new KalmanEstimatorBuilder().
addPropagationConfiguration(propagatorBuilder, new ConstantProcessNoise(initialP, Q)).
build();
// Filter the measurements and check the results
final double expectedDeltaPos = 0.;
final double posEps = 2.33e-8;
final double expectedDeltaVel = 0.;
final double velEps = 6.59e-11;
final double[] expectedsigmasPos = {0.998894, 0.933798, 0.997346};
final double sigmaPosEps = 1e-6;
final double[] expectedSigmasVel = {9.475825e-4, 9.903811e-4, 5.061704e-4};
final double sigmaVelEps = 1e-10;
BrouwerLyddaneEstimationTestUtils.checkKalmanFit(context, kalman, measurements,
refOrbit, positionAngle,
expectedDeltaPos, posEps,
expectedDeltaVel, velEps,
expectedsigmasPos, sigmaPosEps,
expectedSigmasVel, sigmaVelEps);
}
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment