1 package org.orekit.propagation.analytical.tle;
2
3 import org.hipparchus.linear.RealMatrix;
4 import org.hipparchus.util.FastMath;
5 import org.junit.jupiter.api.Assertions;
6 import org.junit.jupiter.api.BeforeEach;
7 import org.junit.jupiter.api.Test;
8 import org.orekit.Utils;
9 import org.orekit.orbits.OrbitType;
10 import org.orekit.orbits.PositionAngleType;
11 import org.orekit.propagation.SpacecraftState;
12 import org.orekit.time.AbsoluteDate;
13 import org.orekit.utils.ParameterDriver;
14 import org.orekit.utils.ParameterDriversList;
15
16 public class TLEParametersDerivativesTest {
17
18
19 private TLE tleSPOT;
20
21 @BeforeEach
22 public void setUp() {
23 Utils.setDataRoot("regular-data");
24
25 String line1SPOT = "1 22823U 93061A 03339.49496229 .00000173 00000-0 10336-3 0 133";
26 String line2SPOT = "2 22823 98.4132 359.2998 0017888 100.4310 259.8872 14.18403464527664";
27 tleSPOT = new TLE(line1SPOT, line2SPOT);
28 }
29
30 @Test
31 public void testBStarEstimation() {
32 doTestParametersDerivatives(TLE.B_STAR, 5.16e-3, tleSPOT);
33 }
34
35 @Test
36 public void testNoEstimatedParameters() {
37
38 TLEPropagator propagator = TLEPropagator.selectExtrapolator(tleSPOT);
39 final SpacecraftState initialState = propagator.getInitialState();
40 final double[] stateVector = new double[6];
41 OrbitType.CARTESIAN.mapOrbitToArray(initialState.getOrbit(), PositionAngleType.MEAN, stateVector, null);
42 TLEHarvester harvester = (TLEHarvester) propagator.setupMatricesComputation("stm", null, null);
43 harvester.freezeColumnsNames();
44 RealMatrix dYdP = harvester.getParametersJacobian(initialState);
45 Assertions.assertNull(dYdP);
46 }
47
48 private void doTestParametersDerivatives(String parameterName, double tolerance, TLE tle) {
49
50
51 ParameterDriversList bound = new ParameterDriversList();
52
53 for (final ParameterDriver driver : tle.getParametersDrivers()) {
54 if (driver.getName().equals(parameterName)) {
55 driver.setSelected(true);
56 bound.add(driver);
57 } else {
58 driver.setSelected(false);
59 }
60 }
61
62 TLEPropagator propagator = TLEPropagator.selectExtrapolator(tle);
63 final SpacecraftState initialState = propagator.getInitialState();
64 final double[] stateVector = new double[6];
65 OrbitType.CARTESIAN.mapOrbitToArray(initialState.getOrbit(), PositionAngleType.MEAN, stateVector, null);
66 final AbsoluteDate target = initialState.getDate().shiftedBy(initialState.getKeplerianPeriod());
67 TLEHarvester harvester = (TLEHarvester) propagator.setupMatricesComputation("stm", null, null);
68 harvester.freezeColumnsNames();
69 RealMatrix dYdP = harvester.getParametersJacobian(initialState);
70 for (int i = 0; i < 6; ++i) {
71 for (int j = 0; j < 1; ++j) {
72 Assertions.assertEquals(0.0, dYdP.getEntry(i, j), tolerance);
73 }
74 }
75 final SpacecraftState finalState = propagator.propagate(target);
76 dYdP = harvester.getParametersJacobian(finalState);
77
78
79 OrbitType orbitType = OrbitType.CARTESIAN;
80 TLEPropagator propagator2 = TLEPropagator.selectExtrapolator(tle);
81 double[][] dYdPRef = new double[6][1];
82
83 ParameterDriver selected = bound.getDrivers().get(0);
84 double p0 = selected.getReferenceValue();
85 double h = selected.getScale();
86 selected.setValue(p0 - 4 * h);
87 propagator2 = TLEPropagator.selectExtrapolator(newTLE(tle, selected.getValue()));
88 SpacecraftState sM4h = propagator2.propagate(target);
89 selected.setValue(p0 - 3 * h);
90 propagator2 = TLEPropagator.selectExtrapolator(newTLE(tle, selected.getValue()));
91 SpacecraftState sM3h = propagator2.propagate(target);
92 selected.setValue(p0 - 2 * h);
93 propagator2 = TLEPropagator.selectExtrapolator(newTLE(tle, selected.getValue()));
94 SpacecraftState sM2h = propagator2.propagate(target);
95 selected.setValue(p0 - 1 * h);
96 propagator2 = TLEPropagator.selectExtrapolator(newTLE(tle, selected.getValue()));
97 SpacecraftState sM1h = propagator2.propagate(target);
98 selected.setValue(p0 + 1 * h);
99 propagator2 = TLEPropagator.selectExtrapolator(newTLE(tle, selected.getValue()));
100 SpacecraftState sP1h = propagator2.propagate(target);
101 selected.setValue(p0 + 2 * h);
102 propagator2 = TLEPropagator.selectExtrapolator(newTLE(tle, selected.getValue()));
103 SpacecraftState sP2h = propagator2.propagate(target);
104 selected.setValue(p0 + 3 * h);
105 propagator2 = TLEPropagator.selectExtrapolator(newTLE(tle, selected.getValue()));
106 SpacecraftState sP3h = propagator2.propagate(target);
107 selected.setValue(p0 + 4 * h);
108 propagator2 = TLEPropagator.selectExtrapolator(newTLE(tle, selected.getValue()));
109 SpacecraftState sP4h = propagator2.propagate(target);
110 fillJacobianColumn(dYdPRef, 0, orbitType, h,
111 sM4h, sM3h, sM2h, sM1h, sP1h, sP2h, sP3h, sP4h);
112
113 for (int i = 0; i < 6; ++i) {
114 Assertions.assertEquals(dYdPRef[i][0], dYdP.getEntry(i, 0), FastMath.abs(tolerance));
115 }
116
117 }
118
119 private void fillJacobianColumn(double[][] jacobian, int column,
120 OrbitType orbitType, double h,
121 SpacecraftState sM4h, SpacecraftState sM3h,
122 SpacecraftState sM2h, SpacecraftState sM1h,
123 SpacecraftState sP1h, SpacecraftState sP2h,
124 SpacecraftState sP3h, SpacecraftState sP4h) {
125 double[] aM4h = stateToArray(sM4h, orbitType)[0];
126 double[] aM3h = stateToArray(sM3h, orbitType)[0];
127 double[] aM2h = stateToArray(sM2h, orbitType)[0];
128 double[] aM1h = stateToArray(sM1h, orbitType)[0];
129 double[] aP1h = stateToArray(sP1h, orbitType)[0];
130 double[] aP2h = stateToArray(sP2h, orbitType)[0];
131 double[] aP3h = stateToArray(sP3h, orbitType)[0];
132 double[] aP4h = stateToArray(sP4h, orbitType)[0];
133 for (int i = 0; i < jacobian.length; ++i) {
134 jacobian[i][column] = ( -3 * (aP4h[i] - aM4h[i]) +
135 32 * (aP3h[i] - aM3h[i]) -
136 168 * (aP2h[i] - aM2h[i]) +
137 672 * (aP1h[i] - aM1h[i])) / (840 * h);
138 }
139 }
140
141 private double[][] stateToArray(SpacecraftState state, OrbitType orbitType) {
142 double[][] array = new double[2][6];
143
144 orbitType.mapOrbitToArray(state.getOrbit(), PositionAngleType.MEAN, array[0], array[1]);
145 return array;
146 }
147
148 private TLE newTLE(final TLE template, final double newBStar) {
149 return new TLE(template.getSatelliteNumber(), template.getClassification(),
150 template.getLaunchYear(), template.getLaunchNumber(), template.getLaunchPiece(),
151 template.getEphemerisType(), template.getElementNumber(), template.getDate(),
152 template.getMeanMotion(), template.getMeanMotionFirstDerivative(), template.getMeanMotionSecondDerivative(),
153 template.getE(), template.getI(), template.getPerigeeArgument(),
154 template.getRaan(), template.getMeanAnomaly(), template.getRevolutionNumberAtEpoch(),
155 newBStar);
156 }
157
158 }