Commit 2fd0ba84 authored by Bryan Cazabonne's avatar Bryan Cazabonne
Browse files

Merge branch 'issue-756' into 'develop'

Added new method signature in IodGooding using AngularRaDec measurement.

Closes #756

See merge request orekit/orekit!136
parents 87436ecf f636ae58
Pipeline #931 failed with stage
......@@ -21,6 +21,9 @@
</properties>
<body>
<release version="11.0" date="TBD" description="TBD">
<action dev="bryan" type="add" issue="756">
Added new method signature in IodGooding using AngularRaDec measurement.
</action>
<action dev="thomas" type="fix" issue="688">
Fixed ignored fields from TLE template in TLEPropagatorBuilder.
</action>
......
......@@ -18,6 +18,7 @@ package org.orekit.estimation.iod;
import org.hipparchus.geometry.euclidean.threed.Vector3D;
import org.hipparchus.util.FastMath;
import org.orekit.estimation.measurements.AngularRaDec;
import org.orekit.frames.Frame;
import org.orekit.orbits.KeplerianOrbit;
import org.orekit.time.AbsoluteDate;
......@@ -121,6 +122,48 @@ public class IodGooding {
return rho3 * R;
}
/** Orbit got from three angular observations.
* @param frame inertial frame for observer coordinates and orbit estimate
* @param raDec1 first angular observation
* @param raDec2 second angular observation
* @param raDec3 third angular observation
* @param rho1init initial guess of the range problem. range 1, in meters
* @param rho3init initial guess of the range problem. range 3, in meters
* @param nRev number of complete revolutions between observation 1 and 3
* @param direction true if posigrade (short way)
* @return an estimate of the Keplerian orbit
* @since 11.0
*/
public KeplerianOrbit estimate(final Frame frame, final AngularRaDec raDec1,
final AngularRaDec raDec2, final AngularRaDec raDec3,
final double rho1init, final double rho3init,
final int nRev, final boolean direction) {
return estimate(frame,
stationPosition(frame, raDec1),
stationPosition(frame, raDec2),
stationPosition(frame, raDec3),
lineOfSight(raDec1), raDec1.getDate(),
lineOfSight(raDec2), raDec2.getDate(),
lineOfSight(raDec3), raDec3.getDate(),
rho1init, rho3init, nRev, direction);
}
/** Orbit got from three angular observations.
* @param frame inertial frame for observer coordinates and orbit estimate
* @param raDec1 first angular observation
* @param raDec2 second angular observation
* @param raDec3 third angular observation
* @param rho1init initial guess of the range problem. range 1, in meters
* @param rho3init initial guess of the range problem. range 3, in meters
* @return an estimate of the Keplerian orbit
* @since 11.0
*/
public KeplerianOrbit estimate(final Frame frame, final AngularRaDec raDec1,
final AngularRaDec raDec2, final AngularRaDec raDec3,
final double rho1init, final double rho3init) {
return estimate(frame, raDec1, raDec2, raDec3, rho1init, rho3init, 0, true);
}
/** Orbit got from Observed Three Lines of Sight (angles only).
* @param frame inertial frame for observer coordinates and orbit estimate
* @param O1 Observer position 1
......@@ -588,4 +631,43 @@ public class IodGooding {
return orbit.shiftedBy(tau).getPVCoordinates().getPosition();
}
/**
* Calculates the line of sight vector.
* @param alpha right ascension angle, in radians
* @param delta declination angle, in radians
* @return the line of sight vector
* @since 11.0
*/
public static Vector3D lineOfSight(final double alpha, final double delta) {
return new Vector3D(FastMath.cos(delta) * FastMath.cos(alpha),
FastMath.cos(delta) * FastMath.sin(alpha),
FastMath.sin(delta));
}
/**
* Calculate the line of sight vector from an AngularRaDec measurement.
* @param raDec measurement
* @return the line of sight vector
* @since 11.0
*/
public static Vector3D lineOfSight(final AngularRaDec raDec) {
// Observed values
final double[] observed = raDec.getObservedValue();
// Return
return lineOfSight(observed[0], observed[1]);
}
/**
* Get the station position from the right ascension / declination measurement.
* @param frame inertial frame for station posiiton and orbit estimate
* @param raDec measurement
* @return the station position
*/
private static Vector3D stationPosition(final Frame frame, final AngularRaDec raDec) {
return raDec.getStation().getBaseFrame().getPVCoordinates(raDec.getDate(), frame).getPosition();
}
}
......@@ -179,7 +179,7 @@ public class IodLaplace {
* @return the line of sight vector
* @since 11.0
*/
private static Vector3D lineOfSight(final AngularRaDec raDec) {
public static Vector3D lineOfSight(final AngularRaDec raDec) {
// Observed values
final double[] observed = raDec.getObservedValue();
......
......@@ -24,6 +24,8 @@ import org.junit.Assert;
import org.junit.Test;
import org.orekit.estimation.Context;
import org.orekit.estimation.EstimationTestUtils;
import org.orekit.estimation.measurements.AngularRaDec;
import org.orekit.estimation.measurements.AngularRaDecMeasurementCreator;
import org.orekit.estimation.measurements.ObservedMeasurement;
import org.orekit.estimation.measurements.PVMeasurementCreator;
import org.orekit.frames.Frame;
......@@ -124,4 +126,62 @@ public class IodGoodingTest {
Assert.assertEquals(13950296.64852, iod.getRange3(), 1.0e-3);
}
@Test
public void testIssue756() {
final Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
final double mu = context.initialOrbit.getMu();
final Frame frame = context.initialOrbit.getFrame();
final NumericalPropagatorBuilder propagatorBuilder =
context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true,
1.0e-6, 60.0, 0.001);
// create perfect range measurements
final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
propagatorBuilder);
final List<ObservedMeasurement<?>> measurements =
EstimationTestUtils.createMeasurements(propagator,
new AngularRaDecMeasurementCreator(context),
0.0, 1.0, 60.0);
// Angular measurements
final AngularRaDec raDec1 = (AngularRaDec) measurements.get(0);
final AngularRaDec raDec2 = (AngularRaDec) measurements.get(20);
final AngularRaDec raDec3 = (AngularRaDec) measurements.get(40);
// Range estimations
final double rhoInit1 = 1.3127847998082995E7;
final double rhoInit3 = 1.3950296648518201E7;
// instantiate the IOD method
final IodGooding iod = new IodGooding(mu);
// TODO convert angular to line of sight and compare both computations
final KeplerianOrbit orbit1 = iod.estimate(frame, raDec1, raDec2, raDec3, rhoInit1, rhoInit3);
final KeplerianOrbit orbit2 = iod.estimate(frame,
stationPosition(frame, raDec1),
stationPosition(frame, raDec2),
stationPosition(frame, raDec3),
IodGooding.lineOfSight(raDec1), raDec1.getDate(),
IodGooding.lineOfSight(raDec2), raDec2.getDate(),
IodGooding.lineOfSight(raDec3), raDec3.getDate(),
rhoInit1, rhoInit3);
Assert.assertEquals(orbit1.getA(), orbit2.getA(), 1.0e-6 * orbit2.getA());
Assert.assertEquals(orbit1.getE(), orbit2.getE(), 1.0e-6 * orbit2.getE());
Assert.assertEquals(orbit1.getI(), orbit2.getI(), 1.0e-6 * orbit2.getI());
Assert.assertEquals(orbit1.getRightAscensionOfAscendingNode(), orbit2.getRightAscensionOfAscendingNode(), 1.0e-6 * orbit2.getRightAscensionOfAscendingNode());
Assert.assertEquals(orbit1.getPerigeeArgument(), orbit2.getPerigeeArgument(), 1.0e-6 * orbit2.getPerigeeArgument());
Assert.assertEquals(orbit1.getMeanAnomaly(), orbit2.getMeanAnomaly(), 1.0e-6 * orbit2.getMeanAnomaly());
}
private static Vector3D stationPosition(final Frame frame, final AngularRaDec raDec) {
return raDec.getStation().getBaseFrame().getPVCoordinates(raDec.getDate(), frame).getPosition();
}
}
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