1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17 package org.orekit.orbits;
18
19 import org.hamcrest.MatcherAssert;
20 import org.hipparchus.analysis.UnivariateFunction;
21 import org.hipparchus.analysis.differentiation.DSFactory;
22 import org.hipparchus.analysis.differentiation.FiniteDifferencesDifferentiator;
23 import org.hipparchus.analysis.differentiation.UnivariateDifferentiableFunction;
24 import org.hipparchus.geometry.euclidean.threed.Vector3D;
25 import org.hipparchus.linear.MatrixUtils;
26 import org.hipparchus.linear.RealMatrixPreservingVisitor;
27 import org.hipparchus.util.FastMath;
28 import org.hipparchus.util.MathUtils;
29 import org.junit.jupiter.api.AfterEach;
30 import org.junit.jupiter.api.Assertions;
31 import org.junit.jupiter.api.BeforeEach;
32 import org.junit.jupiter.api.Test;
33 import org.orekit.Utils;
34 import org.orekit.errors.OrekitException;
35 import org.orekit.errors.OrekitIllegalArgumentException;
36 import org.orekit.errors.OrekitMessages;
37 import org.orekit.frames.Frame;
38 import org.orekit.frames.FramesFactory;
39 import org.orekit.frames.Transform;
40 import org.orekit.propagation.analytical.EcksteinHechlerPropagator;
41 import org.orekit.time.AbsoluteDate;
42 import org.orekit.time.TimeScalesFactory;
43 import org.orekit.utils.Constants;
44 import org.orekit.utils.PVCoordinates;
45 import org.orekit.utils.TimeStampedPVCoordinates;
46
47 import java.io.ByteArrayInputStream;
48 import java.io.ByteArrayOutputStream;
49 import java.io.IOException;
50 import java.io.ObjectInputStream;
51 import java.io.ObjectOutputStream;
52 import java.util.ArrayList;
53 import java.util.List;
54 import java.util.function.Function;
55
56 import static org.orekit.OrekitMatchers.relativelyCloseTo;
57
58 public class KeplerianOrbitTest {
59
60
61 private AbsoluteDate date;
62
63
64 private double mu;
65
66 @Test
67 public void testKeplerianToKeplerian() {
68
69
70 KeplerianOrbit kep =
71 new KeplerianOrbit(24464560.0, 0.7311, 0.122138, 3.10686, 1.00681,
72 0.048363, PositionAngle.MEAN,
73 FramesFactory.getEME2000(), date, mu);
74
75 Vector3D pos = kep.getPVCoordinates().getPosition();
76 Vector3D vit = kep.getPVCoordinates().getVelocity();
77
78 KeplerianOrbit param = new KeplerianOrbit(new PVCoordinates(pos, vit),
79 FramesFactory.getEME2000(), date, mu);
80 Assertions.assertEquals(param.getA(), kep.getA(), Utils.epsilonTest * kep.getA());
81 Assertions.assertEquals(param.getE(), kep.getE(), Utils.epsilonE * FastMath.abs(kep.getE()));
82 Assertions.assertEquals(MathUtils.normalizeAngle(param.getI(), kep.getI()), kep.getI(), Utils.epsilonAngle * FastMath.abs(kep.getI()));
83 Assertions.assertEquals(MathUtils.normalizeAngle(param.getPerigeeArgument(), kep.getPerigeeArgument()), kep.getPerigeeArgument(), Utils.epsilonAngle * FastMath.abs(kep.getPerigeeArgument()));
84 Assertions.assertEquals(MathUtils.normalizeAngle(param.getRightAscensionOfAscendingNode(), kep.getRightAscensionOfAscendingNode()), kep.getRightAscensionOfAscendingNode(), Utils.epsilonAngle * FastMath.abs(kep.getRightAscensionOfAscendingNode()));
85 Assertions.assertEquals(MathUtils.normalizeAngle(param.getMeanAnomaly(), kep.getMeanAnomaly()), kep.getMeanAnomaly(), Utils.epsilonAngle * FastMath.abs(kep.getMeanAnomaly()));
86
87
88 KeplerianOrbit kepCir =
89 new KeplerianOrbit(24464560.0, 0.0, 0.122138, 3.10686, 1.00681,
90 0.048363, PositionAngle.MEAN,
91 FramesFactory.getEME2000(), date, mu);
92
93 Vector3D posCir = kepCir.getPVCoordinates().getPosition();
94 Vector3D vitCir = kepCir.getPVCoordinates().getVelocity();
95
96 KeplerianOrbit paramCir = new KeplerianOrbit(new PVCoordinates(posCir, vitCir),
97 FramesFactory.getEME2000(), date, mu);
98 Assertions.assertEquals(paramCir.getA(), kepCir.getA(), Utils.epsilonTest * kepCir.getA());
99 Assertions.assertEquals(paramCir.getE(), kepCir.getE(), Utils.epsilonE * FastMath.max(1., FastMath.abs(kepCir.getE())));
100 Assertions.assertEquals(MathUtils.normalizeAngle(paramCir.getI(), kepCir.getI()), kepCir.getI(), Utils.epsilonAngle * FastMath.abs(kepCir.getI()));
101 Assertions.assertEquals(MathUtils.normalizeAngle(paramCir.getLM(), kepCir.getLM()), kepCir.getLM(), Utils.epsilonAngle * FastMath.abs(kepCir.getLM()));
102 Assertions.assertEquals(MathUtils.normalizeAngle(paramCir.getLE(), kepCir.getLE()), kepCir.getLE(), Utils.epsilonAngle * FastMath.abs(kepCir.getLE()));
103 Assertions.assertEquals(MathUtils.normalizeAngle(paramCir.getLv(), kepCir.getLv()), kepCir.getLv(), Utils.epsilonAngle * FastMath.abs(kepCir.getLv()));
104
105
106 KeplerianOrbit kepHyp =
107 new KeplerianOrbit(-24464560.0, 1.7311, 0.122138, 3.10686, 1.00681,
108 0.048363, PositionAngle.MEAN,
109 FramesFactory.getEME2000(), date, mu);
110
111 Vector3D posHyp = kepHyp.getPVCoordinates().getPosition();
112 Vector3D vitHyp = kepHyp.getPVCoordinates().getVelocity();
113
114 KeplerianOrbit paramHyp = new KeplerianOrbit(new PVCoordinates(posHyp, vitHyp),
115 FramesFactory.getEME2000(), date, mu);
116 Assertions.assertEquals(paramHyp.getA(), kepHyp.getA(), Utils.epsilonTest * FastMath.abs(kepHyp.getA()));
117 Assertions.assertEquals(paramHyp.getE(), kepHyp.getE(), Utils.epsilonE * FastMath.abs(kepHyp.getE()));
118 Assertions.assertEquals(MathUtils.normalizeAngle(paramHyp.getI(), kepHyp.getI()), kepHyp.getI(), Utils.epsilonAngle * FastMath.abs(kepHyp.getI()));
119 Assertions.assertEquals(MathUtils.normalizeAngle(paramHyp.getPerigeeArgument(), kepHyp.getPerigeeArgument()), kepHyp.getPerigeeArgument(), Utils.epsilonAngle * FastMath.abs(kepHyp.getPerigeeArgument()));
120 Assertions.assertEquals(MathUtils.normalizeAngle(paramHyp.getRightAscensionOfAscendingNode(), kepHyp.getRightAscensionOfAscendingNode()), kepHyp.getRightAscensionOfAscendingNode(), Utils.epsilonAngle * FastMath.abs(kepHyp.getRightAscensionOfAscendingNode()));
121 Assertions.assertEquals(MathUtils.normalizeAngle(paramHyp.getMeanAnomaly(), kepHyp.getMeanAnomaly()), kepHyp.getMeanAnomaly(), Utils.epsilonAngle * FastMath.abs(kepHyp.getMeanAnomaly()));
122
123 }
124
125 @Test
126 public void testKeplerianToCartesian() {
127
128 KeplerianOrbit kep =
129 new KeplerianOrbit(24464560.0, 0.7311, 0.122138, 3.10686, 1.00681,
130 0.048363, PositionAngle.MEAN,
131 FramesFactory.getEME2000(), date, mu);
132
133 Vector3D pos = kep.getPVCoordinates().getPosition();
134 Vector3D vit = kep.getPVCoordinates().getVelocity();
135 Assertions.assertEquals(-0.107622532467967e+07, pos.getX(), Utils.epsilonTest * FastMath.abs(pos.getX()));
136 Assertions.assertEquals(-0.676589636432773e+07, pos.getY(), Utils.epsilonTest * FastMath.abs(pos.getY()));
137 Assertions.assertEquals(-0.332308783350379e+06, pos.getZ(), Utils.epsilonTest * FastMath.abs(pos.getZ()));
138
139 Assertions.assertEquals( 0.935685775154103e+04, vit.getX(), Utils.epsilonTest * FastMath.abs(vit.getX()));
140 Assertions.assertEquals(-0.331234775037644e+04, vit.getY(), Utils.epsilonTest * FastMath.abs(vit.getY()));
141 Assertions.assertEquals(-0.118801577532701e+04, vit.getZ(), Utils.epsilonTest * FastMath.abs(vit.getZ()));
142 }
143
144 @Test
145 public void testKeplerianToEquinoctial() {
146
147 KeplerianOrbit kep =
148 new KeplerianOrbit(24464560.0, 0.7311, 0.122138, 3.10686, 1.00681,
149 0.048363, PositionAngle.MEAN,
150 FramesFactory.getEME2000(), date, mu);
151
152 Assertions.assertEquals(24464560.0, kep.getA(), Utils.epsilonTest * kep.getA());
153 Assertions.assertEquals(-0.412036802887626, kep.getEquinoctialEx(), Utils.epsilonE * FastMath.abs(kep.getE()));
154 Assertions.assertEquals(-0.603931190671706, kep.getEquinoctialEy(), Utils.epsilonE * FastMath.abs(kep.getE()));
155 Assertions.assertEquals(MathUtils.normalizeAngle(2*FastMath.asin(FastMath.sqrt((FastMath.pow(0.652494417368829e-01, 2)+FastMath.pow(0.103158450084864, 2))/4.)), kep.getI()), kep.getI(), Utils.epsilonAngle * FastMath.abs(kep.getI()));
156 Assertions.assertEquals(MathUtils.normalizeAngle(0.416203300000000e+01, kep.getLM()), kep.getLM(), Utils.epsilonAngle * FastMath.abs(kep.getLM()));
157
158 }
159
160 @Test
161 public void testAnomaly() {
162
163 Vector3D position = new Vector3D(7.0e6, 1.0e6, 4.0e6);
164 Vector3D velocity = new Vector3D(-500.0, 8000.0, 1000.0);
165 double mu = 3.9860047e14;
166
167 KeplerianOrbit p = new KeplerianOrbit(new PVCoordinates(position, velocity),
168 FramesFactory.getEME2000(), date, mu);
169
170
171 double e = p.getE();
172 double eRatio = FastMath.sqrt((1 - e) / (1 + e));
173
174 double v = 1.1;
175
176 double E = 2 * FastMath.atan(eRatio * FastMath.tan(v / 2));
177 double M = E - e * FastMath.sin(E);
178
179 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
180 p.getRightAscensionOfAscendingNode(), v , PositionAngle.TRUE,
181 p.getFrame(), p.getDate(), p.getMu());
182 Assertions.assertEquals(p.getTrueAnomaly(), v, Utils.epsilonAngle * FastMath.abs(v));
183 Assertions.assertEquals(p.getEccentricAnomaly(), E, Utils.epsilonAngle * FastMath.abs(E));
184 Assertions.assertEquals(p.getMeanAnomaly(), M, Utils.epsilonAngle * FastMath.abs(M));
185 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
186 p.getRightAscensionOfAscendingNode(), 0 , PositionAngle.TRUE,
187 p.getFrame(), p.getDate(), p.getMu());
188
189 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
190 p.getRightAscensionOfAscendingNode(), E , PositionAngle.ECCENTRIC,
191 p.getFrame(), p.getDate(), p.getMu());
192 Assertions.assertEquals(p.getTrueAnomaly(), v, Utils.epsilonAngle * FastMath.abs(v));
193 Assertions.assertEquals(p.getEccentricAnomaly(), E, Utils.epsilonAngle * FastMath.abs(E));
194 Assertions.assertEquals(p.getMeanAnomaly(), M, Utils.epsilonAngle * FastMath.abs(M));
195 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
196 p.getRightAscensionOfAscendingNode(), 0 , PositionAngle.TRUE,
197 p.getFrame(), p.getDate(), p.getMu());
198
199 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
200 p.getRightAscensionOfAscendingNode(), M, PositionAngle.MEAN,
201 p.getFrame(), p.getDate(), p.getMu());
202 Assertions.assertEquals(p.getTrueAnomaly(), v, Utils.epsilonAngle * FastMath.abs(v));
203 Assertions.assertEquals(p.getEccentricAnomaly(), E, Utils.epsilonAngle * FastMath.abs(E));
204 Assertions.assertEquals(p.getMeanAnomaly(), M, Utils.epsilonAngle * FastMath.abs(M));
205
206
207 p = new KeplerianOrbit(p.getA(), 0, p.getI(), p.getPerigeeArgument(),
208 p.getRightAscensionOfAscendingNode(), p.getLv() , PositionAngle.TRUE,
209 p.getFrame(), p.getDate(), p.getMu());
210
211 E = v;
212 M = E;
213
214 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
215 p.getRightAscensionOfAscendingNode(), v , PositionAngle.TRUE,
216 p.getFrame(), p.getDate(), p.getMu());
217 Assertions.assertEquals(p.getTrueAnomaly(), v, Utils.epsilonAngle * FastMath.abs(v));
218 Assertions.assertEquals(p.getEccentricAnomaly(), E, Utils.epsilonAngle * FastMath.abs(E));
219 Assertions.assertEquals(p.getMeanAnomaly(), M, Utils.epsilonAngle * FastMath.abs(M));
220 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
221 p.getRightAscensionOfAscendingNode(), 0 , PositionAngle.TRUE,
222 p.getFrame(), p.getDate(), p.getMu());
223
224 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
225 p.getRightAscensionOfAscendingNode(), E , PositionAngle.ECCENTRIC, p.getFrame(), p.getDate(), p.getMu());
226 Assertions.assertEquals(p.getTrueAnomaly(), v, Utils.epsilonAngle * FastMath.abs(v));
227 Assertions.assertEquals(p.getEccentricAnomaly(), E, Utils.epsilonAngle * FastMath.abs(E));
228 Assertions.assertEquals(p.getMeanAnomaly(), M, Utils.epsilonAngle * FastMath.abs(M));
229 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
230 p.getRightAscensionOfAscendingNode(), 0 , PositionAngle.TRUE,
231 p.getFrame(), p.getDate(), p.getMu());
232
233 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
234 p.getRightAscensionOfAscendingNode(), M, PositionAngle.MEAN,
235 p.getFrame(), p.getDate(), p.getMu());
236 Assertions.assertEquals(p.getTrueAnomaly(), v, Utils.epsilonAngle * FastMath.abs(v));
237 Assertions.assertEquals(p.getEccentricAnomaly(), E, Utils.epsilonAngle * FastMath.abs(E));
238 Assertions.assertEquals(p.getMeanAnomaly(), M, Utils.epsilonAngle * FastMath.abs(M));
239
240 }
241
242 @Test
243 public void testPositionVelocityNorms() {
244 double mu = 3.9860047e14;
245
246
247 KeplerianOrbit p =
248 new KeplerianOrbit(24464560.0, 0.7311, 2.1, 3.10686, 1.00681,
249 0.67, PositionAngle.TRUE,
250 FramesFactory.getEME2000(), date, mu);
251
252 double e = p.getE();
253 double v = p.getTrueAnomaly();
254 double ksi = 1 + e * FastMath.cos(v);
255 double nu = e * FastMath.sin(v);
256 double epsilon = FastMath.sqrt((1 - e) * (1 + e));
257
258 double a = p.getA();
259 double na = FastMath.sqrt(mu / a);
260
261
262 Assertions.assertEquals(a * epsilon * epsilon / ksi,
263 p.getPVCoordinates().getPosition().getNorm(),
264 Utils.epsilonTest * FastMath.abs(p.getPVCoordinates().getPosition().getNorm()));
265
266
267 Assertions.assertEquals(na * FastMath.sqrt(ksi * ksi + nu * nu) / epsilon,
268 p.getPVCoordinates().getVelocity().getNorm(),
269 Utils.epsilonTest * FastMath.abs(p.getPVCoordinates().getVelocity().getNorm()));
270
271
272
273 KeplerianOrbit pCirEqua =
274 new KeplerianOrbit(24464560.0, 0.1e-10, 0.1e-8, 3.10686, 1.00681,
275 0.67, PositionAngle.TRUE,
276 FramesFactory.getEME2000(), date, mu);
277
278 e = pCirEqua.getE();
279 v = pCirEqua.getTrueAnomaly();
280 ksi = 1 + e * FastMath.cos(v);
281 nu = e * FastMath.sin(v);
282 epsilon = FastMath.sqrt((1 - e) * (1 + e));
283
284 a = pCirEqua.getA();
285 na = FastMath.sqrt(mu / a);
286
287
288 Assertions.assertEquals(a * epsilon * epsilon / ksi,
289 pCirEqua.getPVCoordinates().getPosition().getNorm(),
290 Utils.epsilonTest * FastMath.abs(pCirEqua.getPVCoordinates().getPosition().getNorm()));
291
292
293 Assertions.assertEquals(na * FastMath.sqrt(ksi * ksi + nu * nu) / epsilon,
294 pCirEqua.getPVCoordinates().getVelocity().getNorm(),
295 Utils.epsilonTest * FastMath.abs(pCirEqua.getPVCoordinates().getVelocity().getNorm()));
296 }
297
298 @Test
299 public void testGeometry() {
300 double mu = 3.9860047e14;
301
302
303 KeplerianOrbit p =
304 new KeplerianOrbit(24464560.0, 0.7311, 2.1, 3.10686, 1.00681,
305 0.67, PositionAngle.TRUE,
306 FramesFactory.getEME2000(), date, mu);
307
308 Vector3D position = p.getPVCoordinates().getPosition();
309 Vector3D velocity = p.getPVCoordinates().getVelocity();
310 Vector3D momentum = p.getPVCoordinates().getMomentum().normalize();
311
312 double apogeeRadius = p.getA() * (1 + p.getE());
313 double perigeeRadius = p.getA() * (1 - p.getE());
314
315 for (double lv = 0; lv <= 2 * FastMath.PI; lv += 2 * FastMath.PI/100.) {
316 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
317 p.getRightAscensionOfAscendingNode(), lv , PositionAngle.TRUE,
318 p.getFrame(), p.getDate(), p.getMu());
319 position = p.getPVCoordinates().getPosition();
320
321
322 Assertions.assertTrue((position.getNorm() - apogeeRadius) <= ( apogeeRadius * Utils.epsilonTest));
323 Assertions.assertTrue((position.getNorm() - perigeeRadius) >= (- perigeeRadius * Utils.epsilonTest));
324
325 position = position.normalize();
326 velocity = p.getPVCoordinates().getVelocity();
327 velocity = velocity.normalize();
328
329
330
331
332 Assertions.assertTrue(FastMath.abs(Vector3D.dotProduct(position, momentum)) < Utils.epsilonTest);
333
334 Assertions.assertTrue(FastMath.abs(Vector3D.dotProduct(velocity, momentum)) < Utils.epsilonTest);
335
336 }
337
338
339 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
340 p.getRightAscensionOfAscendingNode(), 0 , PositionAngle.TRUE, p.getFrame(), p.getDate(), p.getMu());
341 Assertions.assertEquals(p.getPVCoordinates().getPosition().getNorm(), perigeeRadius, perigeeRadius * Utils.epsilonTest);
342
343 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
344 p.getRightAscensionOfAscendingNode(), FastMath.PI , PositionAngle.TRUE, p.getFrame(), p.getDate(), p.getMu());
345 Assertions.assertEquals(p.getPVCoordinates().getPosition().getNorm(), apogeeRadius, apogeeRadius * Utils.epsilonTest);
346
347
348
349 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
350 p.getRightAscensionOfAscendingNode(), FastMath.PI - p.getPerigeeArgument() , PositionAngle.TRUE,
351 p.getFrame(), p.getDate(), p.getMu());
352 Assertions.assertTrue(FastMath.abs(p.getPVCoordinates().getPosition().getZ()) < p.getPVCoordinates().getPosition().getNorm() * Utils.epsilonTest);
353 Assertions.assertTrue(p.getPVCoordinates().getVelocity().getZ() < 0);
354
355
356 p = new KeplerianOrbit(p.getA(), p.getE(), p.getI(), p.getPerigeeArgument(),
357 p.getRightAscensionOfAscendingNode(), 2.0 * FastMath.PI - p.getPerigeeArgument() , PositionAngle.TRUE,
358 p.getFrame(), p.getDate(), p.getMu());
359 Assertions.assertTrue(FastMath.abs(p.getPVCoordinates().getPosition().getZ()) < p.getPVCoordinates().getPosition().getNorm() * Utils.epsilonTest);
360 Assertions.assertTrue(p.getPVCoordinates().getVelocity().getZ() > 0);
361
362
363
364 KeplerianOrbit pCirEqua =
365 new KeplerianOrbit(24464560.0, 0.1e-10, 0.1e-8, 3.10686, 1.00681,
366 0.67, PositionAngle.TRUE, FramesFactory.getEME2000(), date, mu);
367
368 position = pCirEqua.getPVCoordinates().getPosition();
369 velocity = pCirEqua.getPVCoordinates().getVelocity();
370 momentum = Vector3D.crossProduct(position, velocity).normalize();
371
372 apogeeRadius = pCirEqua.getA() * (1 + pCirEqua.getE());
373 perigeeRadius = pCirEqua.getA() * (1 - pCirEqua.getE());
374
375 Assertions.assertEquals(perigeeRadius, apogeeRadius, 1.e+4 * Utils.epsilonTest * apogeeRadius);
376
377 for (double lv = 0; lv <= 2 * FastMath.PI; lv += 2 * FastMath.PI/100.) {
378 pCirEqua = new KeplerianOrbit(pCirEqua.getA(), pCirEqua.getE(), pCirEqua.getI(), pCirEqua.getPerigeeArgument(),
379 pCirEqua.getRightAscensionOfAscendingNode(), lv, PositionAngle.TRUE,
380 pCirEqua.getFrame(), pCirEqua.getDate(), pCirEqua.getMu());
381 position = pCirEqua.getPVCoordinates().getPosition();
382
383
384
385 Assertions.assertTrue((position.getNorm() - apogeeRadius) <= ( apogeeRadius * Utils.epsilonTest));
386 Assertions.assertTrue((position.getNorm() - perigeeRadius) >= (- perigeeRadius * Utils.epsilonTest));
387
388 position = position.normalize();
389 velocity = pCirEqua.getPVCoordinates().getVelocity();
390 velocity = velocity.normalize();
391
392
393
394
395 Assertions.assertTrue(FastMath.abs(Vector3D.dotProduct(position, momentum)) < Utils.epsilonTest);
396
397 Assertions.assertTrue(FastMath.abs(Vector3D.dotProduct(velocity, momentum)) < Utils.epsilonTest);
398
399 }
400 }
401
402 @Test
403 public void testSymmetry() {
404
405
406 Vector3D position = new Vector3D(-4947831., -3765382., -3708221.);
407 Vector3D velocity = new Vector3D(-2079., 5291., -7842.);
408 double mu = 3.9860047e14;
409
410 KeplerianOrbit p = new KeplerianOrbit(new PVCoordinates(position, velocity),
411 FramesFactory.getEME2000(), date, mu);
412 Vector3D positionOffset = p.getPVCoordinates().getPosition().subtract(position);
413 Vector3D velocityOffset = p.getPVCoordinates().getVelocity().subtract(velocity);
414
415 Assertions.assertTrue(positionOffset.getNorm() < Utils.epsilonTest);
416 Assertions.assertTrue(velocityOffset.getNorm() < Utils.epsilonTest);
417
418
419 position = new Vector3D(1742382., -2.440243e7, -0.014517);
420 velocity = new Vector3D(4026.2, 287.479, -3.e-6);
421
422
423 p = new KeplerianOrbit(new PVCoordinates(position, velocity),
424 FramesFactory.getEME2000(), date, mu);
425 positionOffset = p.getPVCoordinates().getPosition().subtract(position);
426 velocityOffset = p.getPVCoordinates().getVelocity().subtract(velocity);
427
428 Assertions.assertTrue(positionOffset.getNorm() < Utils.epsilonTest);
429 Assertions.assertTrue(velocityOffset.getNorm() < Utils.epsilonTest);
430
431 }
432
433 @Test
434 public void testNonInertialFrame() throws IllegalArgumentException {
435 Assertions.assertThrows(IllegalArgumentException.class, () -> {
436 Vector3D position = new Vector3D(-4947831., -3765382., -3708221.);
437 Vector3D velocity = new Vector3D(-2079., 5291., -7842.);
438 PVCoordinates pvCoordinates = new PVCoordinates( position, velocity);
439 new KeplerianOrbit(pvCoordinates,
440 new Frame(FramesFactory.getEME2000(), Transform.IDENTITY, "non-inertial", false),
441 date, mu);
442 });
443 }
444
445 @Test
446 public void testPeriod() {
447 KeplerianOrbit orbit = new KeplerianOrbit(7654321.0, 0.1, 0.2, 0, 0, 0,
448 PositionAngle.TRUE,
449 FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH,
450 mu);
451 Assertions.assertEquals(6664.5521723383589487, orbit.getKeplerianPeriod(), 1.0e-12);
452 Assertions.assertEquals(0.00094277682051291315229, orbit.getKeplerianMeanMotion(), 1.0e-16);
453 }
454
455 @Test
456 public void testHyperbola1() {
457 KeplerianOrbit orbit = new KeplerianOrbit(-10000000.0, 2.5, 0.3, 0, 0, 0.0,
458 PositionAngle.TRUE,
459 FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH,
460 mu);
461 Vector3D perigeeP = orbit.getPVCoordinates().getPosition();
462 Vector3D u = perigeeP.normalize();
463 Vector3D focus1 = Vector3D.ZERO;
464 Vector3D focus2 = new Vector3D(-2 * orbit.getA() * orbit.getE(), u);
465 for (double dt = -5000; dt < 5000; dt += 60) {
466 PVCoordinates pv = orbit.shiftedBy(dt).getPVCoordinates();
467 double d1 = Vector3D.distance(pv.getPosition(), focus1);
468 double d2 = Vector3D.distance(pv.getPosition(), focus2);
469 Assertions.assertEquals(-2 * orbit.getA(), FastMath.abs(d1 - d2), 1.0e-6);
470 KeplerianOrbit rebuilt =
471 new KeplerianOrbit(pv, orbit.getFrame(), orbit.getDate().shiftedBy(dt), mu);
472 Assertions.assertEquals(-10000000.0, rebuilt.getA(), 1.0e-6);
473 Assertions.assertEquals(2.5, rebuilt.getE(), 1.0e-13);
474 }
475 }
476
477 @Test
478 public void testHyperbola2() {
479 KeplerianOrbit orbit = new KeplerianOrbit(-10000000.0, 1.2, 0.3, 0, 0, -1.75,
480 PositionAngle.MEAN,
481 FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH,
482 mu);
483 Vector3D perigeeP = new KeplerianOrbit(orbit.getA(), orbit.getE(), orbit.getI(),
484 orbit.getPerigeeArgument(), orbit.getRightAscensionOfAscendingNode(),
485 0.0, PositionAngle.TRUE, orbit.getFrame(),
486 orbit.getDate(), orbit.getMu()).getPVCoordinates().getPosition();
487 Vector3D u = perigeeP.normalize();
488 Vector3D focus1 = Vector3D.ZERO;
489 Vector3D focus2 = new Vector3D(-2 * orbit.getA() * orbit.getE(), u);
490 for (double dt = -5000; dt < 5000; dt += 60) {
491 PVCoordinates pv = orbit.shiftedBy(dt).getPVCoordinates();
492 double d1 = Vector3D.distance(pv.getPosition(), focus1);
493 double d2 = Vector3D.distance(pv.getPosition(), focus2);
494 Assertions.assertEquals(-2 * orbit.getA(), FastMath.abs(d1 - d2), 1.0e-6);
495 KeplerianOrbit rebuilt =
496 new KeplerianOrbit(pv, orbit.getFrame(), orbit.getDate().shiftedBy(dt), mu);
497 Assertions.assertEquals(-10000000.0, rebuilt.getA(), 1.0e-6);
498 Assertions.assertEquals(1.2, rebuilt.getE(), 1.0e-13);
499 }
500 }
501
502 @Test
503 public void testInconsistentHyperbola() {
504 try {
505 new KeplerianOrbit(+10000000.0, 2.5, 0.3, 0, 0, 0.0,
506 PositionAngle.TRUE,
507 FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH,
508 mu);
509 Assertions.fail("an exception should have been thrown");
510 } catch (OrekitIllegalArgumentException oe) {
511 Assertions.assertEquals(OrekitMessages.ORBIT_A_E_MISMATCH_WITH_CONIC_TYPE, oe.getSpecifier());
512 Assertions.assertEquals(+10000000.0, ((Double) oe.getParts()[0]).doubleValue(), 1.0e-3);
513 Assertions.assertEquals(2.5, ((Double) oe.getParts()[1]).doubleValue(), 1.0e-15);
514 }
515 }
516
517 @Test
518 public void testVeryLargeEccentricity() {
519
520 final Frame eme2000 = FramesFactory.getEME2000();
521 final double meanAnomaly = 1.;
522 final KeplerianOrbit orb0 = new KeplerianOrbit(42600e3, 0.9, 0.00001, 0, 0,
523 FastMath.toRadians(meanAnomaly),
524 PositionAngle.MEAN, eme2000, date, mu);
525
526
527 final Vector3D deltaV = new Vector3D(0.0, 110000.0, 0.0);
528 final PVCoordinates pv1 = new PVCoordinates(orb0.getPVCoordinates().getPosition(),
529 orb0.getPVCoordinates().getVelocity().add(deltaV));
530 final KeplerianOrbit orb1 = new KeplerianOrbit(pv1, eme2000, date, mu);
531
532
533
534 final PVCoordinates pvTested = orb1.shiftedBy(0).getPVCoordinates();
535 final Vector3D pTested = pvTested.getPosition();
536 final Vector3D vTested = pvTested.getVelocity();
537
538 final PVCoordinates pvReference = orb1.getPVCoordinates();
539 final Vector3D pReference = pvReference.getPosition();
540 final Vector3D vReference = pvReference.getVelocity();
541
542 final double threshold = 1.e-15;
543 Assertions.assertEquals(0, pTested.subtract(pReference).getNorm(), threshold * pReference.getNorm());
544 Assertions.assertEquals(0, vTested.subtract(vReference).getNorm(), threshold * vReference.getNorm());
545
546 }
547
548 @Test
549 public void testKeplerEquation() {
550
551 for (double M = -6 * FastMath.PI; M < 6 * FastMath.PI; M += 0.01) {
552 KeplerianOrbit pElliptic =
553 new KeplerianOrbit(24464560.0, 0.7311, 2.1, 3.10686, 1.00681,
554 M, PositionAngle.MEAN,
555 FramesFactory.getEME2000(), date, mu);
556 double E = pElliptic.getEccentricAnomaly();
557 double e = pElliptic.getE();
558 Assertions.assertEquals(M, E - e * FastMath.sin(E), 2.0e-14);
559 }
560
561 for (double M = -6 * FastMath.PI; M < 6 * FastMath.PI; M += 0.01) {
562 KeplerianOrbit pAlmostParabolic =
563 new KeplerianOrbit(24464560.0, 0.9999, 2.1, 3.10686, 1.00681,
564 M, PositionAngle.MEAN,
565 FramesFactory.getEME2000(), date, mu);
566 double E = pAlmostParabolic.getEccentricAnomaly();
567 double e = pAlmostParabolic.getE();
568 Assertions.assertEquals(M, E - e * FastMath.sin(E), 4.0e-13);
569 }
570
571 }
572
573 @Test
574 public void testOutOfRangeV() {
575 Assertions.assertThrows(IllegalArgumentException.class, () -> {
576 new KeplerianOrbit(-7000434.460140012, 1.1999785407363386, 1.3962787004479158,
577 1.3962320168955138, 0.3490728321331678, -2.55593407037698,
578 PositionAngle.TRUE, FramesFactory.getEME2000(),
579 new AbsoluteDate("2000-01-01T12:00:00.391", TimeScalesFactory.getUTC()),
580 3.986004415E14);
581 });
582 }
583
584 @Test
585 public void testNumericalIssue25() {
586 Vector3D position = new Vector3D(3782116.14107698, 416663.11924914, 5875541.62103057);
587 Vector3D velocity = new Vector3D(-6349.7848910501, 288.4061811651, 4066.9366759691);
588 KeplerianOrbit orbit = new KeplerianOrbit(new PVCoordinates(position, velocity),
589 FramesFactory.getEME2000(),
590 new AbsoluteDate("2004-01-01T23:00:00.000",
591 TimeScalesFactory.getUTC()),
592 3.986004415E14);
593 Assertions.assertEquals(0.0, orbit.getE(), 2.0e-14);
594 }
595
596 @Test
597 public void testPerfectlyEquatorial() {
598 Vector3D position = new Vector3D(6957904.3624652653594, 766529.11411558074507, 0);
599 Vector3D velocity = new Vector3D(-7538.2817012412102845, 342.38751001881413381, 0.);
600 KeplerianOrbit orbit = new KeplerianOrbit(new PVCoordinates(position, velocity),
601 FramesFactory.getEME2000(),
602 new AbsoluteDate("2004-01-01T23:00:00.000",
603 TimeScalesFactory.getUTC()),
604 3.986004415E14);
605 Assertions.assertEquals(0.0, orbit.getI(), 2.0e-14);
606 Assertions.assertEquals(0.0, orbit.getRightAscensionOfAscendingNode(), 2.0e-14);
607 }
608
609 @Test
610 public void testJacobianReferenceEllipse() {
611
612 AbsoluteDate dateTca = new AbsoluteDate(2000, 04, 01, 0, 0, 0.000, TimeScalesFactory.getUTC());
613 double mu = 3.986004415e+14;
614 KeplerianOrbit orbKep = new KeplerianOrbit(7000000.0, 0.01, FastMath.toRadians(80.), FastMath.toRadians(80.), FastMath.toRadians(20.),
615 FastMath.toRadians(40.), PositionAngle.MEAN,
616 FramesFactory.getEME2000(), dateTca, mu);
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662 Vector3D pRef = new Vector3D(-3691555.569874833337963, -240330.253992714860942, 5879700.285850423388183);
663 Vector3D vRef = new Vector3D(-5936.229884450408463, -2871.067660163344044, -3786.209549192726627);
664 double[][] jRef = {
665 { -1.0792090588217809, -7.02594292049818631E-002, 1.7189029642216496, -1459.4829009393857, -705.88138246206040, -930.87838644776593 },
666 { -1.31195762636625214E-007, -3.90087231593959271E-008, 4.65917592901869866E-008, -2.02467187867647177E-004, -7.89767994436215424E-005, -2.81639203329454407E-005 },
667 { 4.18334478744371316E-008, -1.14936453412947957E-007, 2.15670500707930151E-008, -2.26450325965329431E-005, 6.22167157217876380E-005, -1.16745469637130306E-005 },
668 { 3.52735168061691945E-006, 3.82555734454450974E-006, 1.34715077236557634E-005, -8.06586262922115264E-003, -6.13725651685311825E-003, -1.71765290503914092E-002 },
669 { 2.48948022169790885E-008, -6.83979069529389238E-008, 1.28344057971888544E-008, 3.86597661353874888E-005, -1.06216834498373629E-004, 1.99308724078785540E-005 },
670 { -3.41911705254704525E-006, -3.75913623359912437E-006, -1.34013845492518465E-005, 8.19851888816422458E-003, 6.16449264680494959E-003, 1.69495878276556648E-002 }
671 };
672
673 PVCoordinates pv = orbKep.getPVCoordinates();
674 Assertions.assertEquals(0, pv.getPosition().subtract(pRef).getNorm(), 1.0e-15 * pRef.getNorm());
675 Assertions.assertEquals(0, pv.getVelocity().subtract(vRef).getNorm(), 1.0e-16 * vRef.getNorm());
676
677 double[][] jacobian = new double[6][6];
678 orbKep.getJacobianWrtCartesian(PositionAngle.MEAN, jacobian);
679
680 for (int i = 0; i < jacobian.length; i++) {
681 double[] row = jacobian[i];
682 double[] rowRef = jRef[i];
683 for (int j = 0; j < row.length; j++) {
684 Assertions.assertEquals(0, (row[j] - rowRef[j]) / rowRef[j], 2.0e-12);
685 }
686 }
687
688 }
689
690 @Test
691 public void testJacobianFinitedifferencesEllipse() {
692
693 AbsoluteDate dateTca = new AbsoluteDate(2000, 04, 01, 0, 0, 0.000, TimeScalesFactory.getUTC());
694 double mu = 3.986004415e+14;
695 KeplerianOrbit orbKep = new KeplerianOrbit(7000000.0, 0.01, FastMath.toRadians(80.), FastMath.toRadians(80.), FastMath.toRadians(20.),
696 FastMath.toRadians(40.), PositionAngle.MEAN,
697 FramesFactory.getEME2000(), dateTca, mu);
698
699 for (PositionAngle type : PositionAngle.values()) {
700 double hP = 2.0;
701 double[][] finiteDiffJacobian = finiteDifferencesJacobian(type, orbKep, hP);
702 double[][] jacobian = new double[6][6];
703 orbKep.getJacobianWrtCartesian(type, jacobian);
704
705 for (int i = 0; i < jacobian.length; i++) {
706 double[] row = jacobian[i];
707 double[] rowRef = finiteDiffJacobian[i];
708 for (int j = 0; j < row.length; j++) {
709 Assertions.assertEquals(0, (row[j] - rowRef[j]) / rowRef[j], 2.0e-7);
710 }
711 }
712
713 double[][] invJacobian = new double[6][6];
714 orbKep.getJacobianWrtParameters(type, invJacobian);
715 MatrixUtils.createRealMatrix(jacobian).
716 multiply(MatrixUtils.createRealMatrix(invJacobian)).
717 walkInRowOrder(new RealMatrixPreservingVisitor() {
718 public void start(int rows, int columns,
719 int startRow, int endRow, int startColumn, int endColumn) {
720 }
721
722 public void visit(int row, int column, double value) {
723 Assertions.assertEquals(row == column ? 1.0 : 0.0, value, 5.0e-9);
724 }
725
726 public double end() {
727 return Double.NaN;
728 }
729 });
730
731 }
732
733 }
734
735 @Test
736 public void testJacobianReferenceHyperbola() {
737
738 AbsoluteDate dateTca = new AbsoluteDate(2000, 04, 01, 0, 0, 0.000, TimeScalesFactory.getUTC());
739 double mu = 3.986004415e+14;
740 KeplerianOrbit orbKep = new KeplerianOrbit(-7000000.0, 1.2, FastMath.toRadians(80.), FastMath.toRadians(80.), FastMath.toRadians(20.),
741 FastMath.toRadians(40.), PositionAngle.MEAN,
742 FramesFactory.getEME2000(), dateTca, mu);
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791 Vector3D pRef = new Vector3D(-7654711.206549182534218, -3460171.872979687992483, -3592374.514463655184954);
792 Vector3D vRef = new Vector3D( -7886.368091820805603, -4359.739012331759113, -7937.060044548694350);
793 double[][] jRef = {
794 { -0.98364725131848019, -0.44463970750901238, -0.46162803814668391, -1938.9443476028839, -1071.8864775981751, -1951.4074832397598 },
795 { -1.10548813242982574E-007, -2.52906747183730431E-008, 7.96500937398593591E-008, -9.70479823470940108E-006, -2.93209076428001017E-005, -1.37434463892791042E-004 },
796 { 8.55737680891616672E-008, -2.35111995522618220E-007, 4.41171797903162743E-008, -8.05235180390949802E-005, 2.21236547547460423E-004, -4.15135455876865407E-005 },
797 { -1.52641427784095578E-007, 1.10250447958827901E-008, 1.21265251605359894E-007, 7.63347077200903542E-005, -3.54738331412232378E-005, -2.31400737283033359E-004 },
798 { 7.86711766048035274E-008, -2.16147281283624453E-007, 4.05585791077187359E-008, -3.56071805267582894E-005, 9.78299244677127374E-005, -1.83571253224293247E-005 },
799 { -2.41488884881911384E-007, -1.00119615610276537E-007, -6.51494225096757969E-008, -2.43295075073248163E-004, -1.43273725071890463E-004, -2.91625510452094873E-004 }
800 };
801
802 PVCoordinates pv = orbKep.getPVCoordinates();
803 Assertions.assertEquals(0, pv.getPosition().subtract(pRef).getNorm() / pRef.getNorm(), 1.0e-16);
804
805 Assertions.assertEquals(0, pv.getVelocity().subtract(vRef).getNorm() / vRef.getNorm(), 3.0e-16);
806
807 double[][] jacobian = new double[6][6];
808 orbKep.getJacobianWrtCartesian(PositionAngle.MEAN, jacobian);
809
810 for (int i = 0; i < jacobian.length; i++) {
811 double[] row = jacobian[i];
812 double[] rowRef = jRef[i];
813 for (int j = 0; j < row.length; j++) {
814 Assertions.assertEquals(0, (row[j] - rowRef[j]) / rowRef[j], 1.0e-14);
815 }
816 }
817
818 }
819
820 @Test
821 public void testJacobianFinitedifferencesHyperbola() {
822
823 AbsoluteDate dateTca = new AbsoluteDate(2000, 04, 01, 0, 0, 0.000, TimeScalesFactory.getUTC());
824 double mu = 3.986004415e+14;
825 KeplerianOrbit orbKep = new KeplerianOrbit(-7000000.0, 1.2, FastMath.toRadians(80.), FastMath.toRadians(80.), FastMath.toRadians(20.),
826 FastMath.toRadians(40.), PositionAngle.MEAN,
827 FramesFactory.getEME2000(), dateTca, mu);
828
829 for (PositionAngle type : PositionAngle.values()) {
830 double hP = 2.0;
831 double[][] finiteDiffJacobian = finiteDifferencesJacobian(type, orbKep, hP);
832 double[][] jacobian = new double[6][6];
833 orbKep.getJacobianWrtCartesian(type, jacobian);
834 for (int i = 0; i < jacobian.length; i++) {
835 double[] row = jacobian[i];
836 double[] rowRef = finiteDiffJacobian[i];
837 for (int j = 0; j < row.length; j++) {
838 Assertions.assertEquals(0, (row[j] - rowRef[j]) / rowRef[j], 3.0e-8);
839 }
840 }
841
842 double[][] invJacobian = new double[6][6];
843 orbKep.getJacobianWrtParameters(type, invJacobian);
844 MatrixUtils.createRealMatrix(jacobian).
845 multiply(MatrixUtils.createRealMatrix(invJacobian)).
846 walkInRowOrder(new RealMatrixPreservingVisitor() {
847 public void start(int rows, int columns,
848 int startRow, int endRow, int startColumn, int endColumn) {
849 }
850
851 public void visit(int row, int column, double value) {
852 Assertions.assertEquals(row == column ? 1.0 : 0.0, value, 2.0e-8);
853 }
854
855 public double end() {
856 return Double.NaN;
857 }
858 });
859
860 }
861
862 }
863
864 private double[][] finiteDifferencesJacobian(PositionAngle type, KeplerianOrbit orbit, double hP)
865 {
866 double[][] jacobian = new double[6][6];
867 for (int i = 0; i < 6; ++i) {
868 fillColumn(type, i, orbit, hP, jacobian);
869 }
870 return jacobian;
871 }
872
873 private void fillColumn(PositionAngle type, int i, KeplerianOrbit orbit, double hP, double[][] jacobian) {
874
875
876
877 Vector3D p = orbit.getPVCoordinates().getPosition();
878 Vector3D v = orbit.getPVCoordinates().getVelocity();
879 double hV = orbit.getMu() * hP / (v.getNorm() * p.getNormSq());
880
881 double h;
882 Vector3D dP = Vector3D.ZERO;
883 Vector3D dV = Vector3D.ZERO;
884 switch (i) {
885 case 0:
886 h = hP;
887 dP = new Vector3D(hP, 0, 0);
888 break;
889 case 1:
890 h = hP;
891 dP = new Vector3D(0, hP, 0);
892 break;
893 case 2:
894 h = hP;
895 dP = new Vector3D(0, 0, hP);
896 break;
897 case 3:
898 h = hV;
899 dV = new Vector3D(hV, 0, 0);
900 break;
901 case 4:
902 h = hV;
903 dV = new Vector3D(0, hV, 0);
904 break;
905 default:
906 h = hV;
907 dV = new Vector3D(0, 0, hV);
908 break;
909 }
910
911 KeplerianOrbit oM4h = new KeplerianOrbit(new PVCoordinates(new Vector3D(1, p, -4, dP), new Vector3D(1, v, -4, dV)),
912 orbit.getFrame(), orbit.getDate(), orbit.getMu());
913 KeplerianOrbit oM3h = new KeplerianOrbit(new PVCoordinates(new Vector3D(1, p, -3, dP), new Vector3D(1, v, -3, dV)),
914 orbit.getFrame(), orbit.getDate(), orbit.getMu());
915 KeplerianOrbit oM2h = new KeplerianOrbit(new PVCoordinates(new Vector3D(1, p, -2, dP), new Vector3D(1, v, -2, dV)),
916 orbit.getFrame(), orbit.getDate(), orbit.getMu());
917 KeplerianOrbit oM1h = new KeplerianOrbit(new PVCoordinates(new Vector3D(1, p, -1, dP), new Vector3D(1, v, -1, dV)),
918 orbit.getFrame(), orbit.getDate(), orbit.getMu());
919 KeplerianOrbit oP1h = new KeplerianOrbit(new PVCoordinates(new Vector3D(1, p, +1, dP), new Vector3D(1, v, +1, dV)),
920 orbit.getFrame(), orbit.getDate(), orbit.getMu());
921 KeplerianOrbit oP2h = new KeplerianOrbit(new PVCoordinates(new Vector3D(1, p, +2, dP), new Vector3D(1, v, +2, dV)),
922 orbit.getFrame(), orbit.getDate(), orbit.getMu());
923 KeplerianOrbit oP3h = new KeplerianOrbit(new PVCoordinates(new Vector3D(1, p, +3, dP), new Vector3D(1, v, +3, dV)),
924 orbit.getFrame(), orbit.getDate(), orbit.getMu());
925 KeplerianOrbit oP4h = new KeplerianOrbit(new PVCoordinates(new Vector3D(1, p, +4, dP), new Vector3D(1, v, +4, dV)),
926 orbit.getFrame(), orbit.getDate(), orbit.getMu());
927
928 jacobian[0][i] = (-3 * (oP4h.getA() - oM4h.getA()) +
929 32 * (oP3h.getA() - oM3h.getA()) -
930 168 * (oP2h.getA() - oM2h.getA()) +
931 672 * (oP1h.getA() - oM1h.getA())) / (840 * h);
932 jacobian[1][i] = (-3 * (oP4h.getE() - oM4h.getE()) +
933 32 * (oP3h.getE() - oM3h.getE()) -
934 168 * (oP2h.getE() - oM2h.getE()) +
935 672 * (oP1h.getE() - oM1h.getE())) / (840 * h);
936 jacobian[2][i] = (-3 * (oP4h.getI() - oM4h.getI()) +
937 32 * (oP3h.getI() - oM3h.getI()) -
938 168 * (oP2h.getI() - oM2h.getI()) +
939 672 * (oP1h.getI() - oM1h.getI())) / (840 * h);
940 jacobian[3][i] = (-3 * (oP4h.getPerigeeArgument() - oM4h.getPerigeeArgument()) +
941 32 * (oP3h.getPerigeeArgument() - oM3h.getPerigeeArgument()) -
942 168 * (oP2h.getPerigeeArgument() - oM2h.getPerigeeArgument()) +
943 672 * (oP1h.getPerigeeArgument() - oM1h.getPerigeeArgument())) / (840 * h);
944 jacobian[4][i] = (-3 * (oP4h.getRightAscensionOfAscendingNode() - oM4h.getRightAscensionOfAscendingNode()) +
945 32 * (oP3h.getRightAscensionOfAscendingNode() - oM3h.getRightAscensionOfAscendingNode()) -
946 168 * (oP2h.getRightAscensionOfAscendingNode() - oM2h.getRightAscensionOfAscendingNode()) +
947 672 * (oP1h.getRightAscensionOfAscendingNode() - oM1h.getRightAscensionOfAscendingNode())) / (840 * h);
948 jacobian[5][i] = (-3 * (oP4h.getAnomaly(type) - oM4h.getAnomaly(type)) +
949 32 * (oP3h.getAnomaly(type) - oM3h.getAnomaly(type)) -
950 168 * (oP2h.getAnomaly(type) - oM2h.getAnomaly(type)) +
951 672 * (oP1h.getAnomaly(type) - oM1h.getAnomaly(type))) / (840 * h);
952
953 }
954
955 @Test
956 public void testInterpolationWithDerivatives() {
957 doTestInterpolation(true,
958 397, 4.01, 4.75e-4, 1.28e-7,
959 2159, 1.05e7, 1.19e-3, 0.773);
960 }
961
962 @Test
963 public void testInterpolationWithoutDerivatives() {
964 doTestInterpolation(false,
965 397, 62.0, 4.75e-4, 2.87e-6,
966 2159, 79365, 1.19e-3, 3.89e-3);
967 }
968
969 private void doTestInterpolation(boolean useDerivatives,
970 double shiftPositionErrorWithin, double interpolationPositionErrorWithin,
971 double shiftEccentricityErrorWithin, double interpolationEccentricityErrorWithin,
972 double shiftPositionErrorSlightlyPast, double interpolationPositionErrorSlightlyPast,
973 double shiftEccentricityErrorSlightlyPast, double interpolationEccentricityErrorSlightlyPast)
974 {
975
976 final double ehMu = 3.9860047e14;
977 final double ae = 6.378137e6;
978 final double c20 = -1.08263e-3;
979 final double c30 = 2.54e-6;
980 final double c40 = 1.62e-6;
981 final double c50 = 2.3e-7;
982 final double c60 = -5.5e-7;
983
984 final AbsoluteDate date = AbsoluteDate.J2000_EPOCH.shiftedBy(584.);
985 final Vector3D position = new Vector3D(3220103., 69623., 6449822.);
986 final Vector3D velocity = new Vector3D(6414.7, -2006., -3180.);
987 final KeplerianOrbit initialOrbit = new KeplerianOrbit(new PVCoordinates(position, velocity),
988 FramesFactory.getEME2000(), date, ehMu);
989
990 EcksteinHechlerPropagator propagator =
991 new EcksteinHechlerPropagator(initialOrbit, ae, ehMu, c20, c30, c40, c50, c60);
992
993
994 List<Orbit> sample = new ArrayList<Orbit>();
995 for (double dt = 0; dt < 300.0; dt += 60.0) {
996 Orbit orbit = propagator.propagate(date.shiftedBy(dt)).getOrbit();
997 if (!useDerivatives) {
998
999 double[] stateVector = new double[6];
1000 orbit.getType().mapOrbitToArray(orbit, PositionAngle.TRUE, stateVector, null);
1001 orbit = orbit.getType().mapArrayToOrbit(stateVector, null, PositionAngle.TRUE,
1002 orbit.getDate(), orbit.getMu(), orbit.getFrame());
1003 }
1004 sample.add(orbit);
1005 }
1006
1007
1008
1009
1010 double maxShiftPositionError = 0;
1011 double maxInterpolationPositionError = 0;
1012 double maxShiftEccentricityError = 0;
1013 double maxInterpolationEccentricityError = 0;
1014 for (double dt = 0; dt < 241.0; dt += 1.0) {
1015 AbsoluteDate t = initialOrbit.getDate().shiftedBy(dt);
1016 Vector3D shiftedP = initialOrbit.shiftedBy(dt).getPVCoordinates().getPosition();
1017 Vector3D interpolatedP = initialOrbit.interpolate(t, sample).getPVCoordinates().getPosition();
1018 Vector3D propagatedP = propagator.propagate(t).getPVCoordinates().getPosition();
1019 double shiftedE = initialOrbit.shiftedBy(dt).getE();
1020 double interpolatedE = initialOrbit.interpolate(t, sample).getE();
1021 double propagatedE = propagator.propagate(t).getE();
1022 maxShiftPositionError = FastMath.max(maxShiftPositionError, shiftedP.subtract(propagatedP).getNorm());
1023 maxInterpolationPositionError = FastMath.max(maxInterpolationPositionError, interpolatedP.subtract(propagatedP).getNorm());
1024 maxShiftEccentricityError = FastMath.max(maxShiftEccentricityError, FastMath.abs(shiftedE - propagatedE));
1025 maxInterpolationEccentricityError = FastMath.max(maxInterpolationEccentricityError, FastMath.abs(interpolatedE - propagatedE));
1026 }
1027 Assertions.assertEquals(shiftPositionErrorWithin, maxShiftPositionError, 0.01 * shiftPositionErrorWithin);
1028 Assertions.assertEquals(interpolationPositionErrorWithin, maxInterpolationPositionError, 0.01 * interpolationPositionErrorWithin);
1029 Assertions.assertEquals(shiftEccentricityErrorWithin, maxShiftEccentricityError, 0.01 * shiftEccentricityErrorWithin);
1030 Assertions.assertEquals(interpolationEccentricityErrorWithin, maxInterpolationEccentricityError, 0.01 * interpolationEccentricityErrorWithin);
1031
1032
1033
1034
1035
1036 maxShiftPositionError = 0;
1037 maxInterpolationPositionError = 0;
1038 maxShiftEccentricityError = 0;
1039 maxInterpolationEccentricityError = 0;
1040 for (double dt = 240; dt < 600; dt += 1.0) {
1041 AbsoluteDate t = initialOrbit.getDate().shiftedBy(dt);
1042 Vector3D shiftedP = initialOrbit.shiftedBy(dt).getPVCoordinates().getPosition();
1043 Vector3D interpolatedP = initialOrbit.interpolate(t, sample).getPVCoordinates().getPosition();
1044 Vector3D propagatedP = propagator.propagate(t).getPVCoordinates().getPosition();
1045 double shiftedE = initialOrbit.shiftedBy(dt).getE();
1046 double interpolatedE = initialOrbit.interpolate(t, sample).getE();
1047 double propagatedE = propagator.propagate(t).getE();
1048 maxShiftPositionError = FastMath.max(maxShiftPositionError, shiftedP.subtract(propagatedP).getNorm());
1049 maxInterpolationPositionError = FastMath.max(maxInterpolationPositionError, interpolatedP.subtract(propagatedP).getNorm());
1050 maxShiftEccentricityError = FastMath.max(maxShiftEccentricityError, FastMath.abs(shiftedE - propagatedE));
1051 maxInterpolationEccentricityError = FastMath.max(maxInterpolationEccentricityError, FastMath.abs(interpolatedE - propagatedE));
1052 }
1053 Assertions.assertEquals(shiftPositionErrorSlightlyPast, maxShiftPositionError, 0.01 * shiftPositionErrorSlightlyPast);
1054 Assertions.assertEquals(interpolationPositionErrorSlightlyPast, maxInterpolationPositionError, 0.01 * interpolationPositionErrorSlightlyPast);
1055 Assertions.assertEquals(shiftEccentricityErrorSlightlyPast, maxShiftEccentricityError, 0.01 * shiftEccentricityErrorSlightlyPast);
1056 Assertions.assertEquals(interpolationEccentricityErrorSlightlyPast, maxInterpolationEccentricityError, 0.01 * interpolationEccentricityErrorSlightlyPast);
1057
1058 }
1059
1060 @Test
1061 public void testPerfectlyEquatorialConversion() {
1062 KeplerianOrbit initial = new KeplerianOrbit(13378000.0, 0.05, 0.0, 0.0, FastMath.PI,
1063 0.0, PositionAngle.MEAN,
1064 FramesFactory.getEME2000(), date,
1065 Constants.EIGEN5C_EARTH_MU);
1066 EquinoctialOrbit equ = (EquinoctialOrbit) OrbitType.EQUINOCTIAL.convertType(initial);
1067 KeplerianOrbit converted = (KeplerianOrbit) OrbitType.KEPLERIAN.convertType(equ);
1068 Assertions.assertEquals(FastMath.PI,
1069 MathUtils.normalizeAngle(converted.getRightAscensionOfAscendingNode() +
1070 converted.getPerigeeArgument(), FastMath.PI),
1071 1.0e-10);
1072 }
1073
1074 @Test
1075 public void testKeplerianDerivatives() {
1076
1077 final KeplerianOrbit orbit = new KeplerianOrbit(new PVCoordinates(new Vector3D(-4947831., -3765382., -3708221.),
1078 new Vector3D(-2079., 5291., -7842.)),
1079 FramesFactory.getEME2000(), date, 3.9860047e14);
1080 final Vector3D p = orbit.getPVCoordinates().getPosition();
1081 final Vector3D v = orbit.getPVCoordinates().getVelocity();
1082 final Vector3D a = orbit.getPVCoordinates().getAcceleration();
1083
1084
1085 Assertions.assertEquals(7.605422, a.getNorm(), 1.0e-6);
1086
1087
1088 Assertions.assertEquals(differentiate(orbit, shifted -> shifted.getPVCoordinates().getPosition().getX()),
1089 orbit.getPVCoordinates().getVelocity().getX(),
1090 3.0e-12 * v.getNorm());
1091 Assertions.assertEquals(differentiate(orbit, shifted -> shifted.getPVCoordinates().getPosition().getY()),
1092 orbit.getPVCoordinates().getVelocity().getY(),
1093 3.0e-12 * v.getNorm());
1094 Assertions.assertEquals(differentiate(orbit, shifted -> shifted.getPVCoordinates().getPosition().getZ()),
1095 orbit.getPVCoordinates().getVelocity().getZ(),
1096 3.0e-12 * v.getNorm());
1097
1098
1099 Assertions.assertEquals(differentiate(orbit, shifted -> shifted.getPVCoordinates().getVelocity().getX()),
1100 orbit.getPVCoordinates().getAcceleration().getX(),
1101 3.0e-12 * a.getNorm());
1102 Assertions.assertEquals(differentiate(orbit, shifted -> shifted.getPVCoordinates().getVelocity().getY()),
1103 orbit.getPVCoordinates().getAcceleration().getY(),
1104 3.0e-12 * a.getNorm());
1105 Assertions.assertEquals(differentiate(orbit, shifted -> shifted.getPVCoordinates().getVelocity().getZ()),
1106 orbit.getPVCoordinates().getAcceleration().getZ(),
1107 3.0e-12 * a.getNorm());
1108
1109
1110 final double r2 = p.getNormSq();
1111 final double r = FastMath.sqrt(r2);
1112 Vector3D keplerianJerk = new Vector3D(-3 * Vector3D.dotProduct(p, v) / r2, a, -a.getNorm() / r, v);
1113 Assertions.assertEquals(differentiate(orbit, shifted -> shifted.getPVCoordinates().getAcceleration().getX()),
1114 keplerianJerk.getX(),
1115 3.0e-12 * keplerianJerk.getNorm());
1116 Assertions.assertEquals(differentiate(orbit, shifted -> shifted.getPVCoordinates().getAcceleration().getY()),
1117 keplerianJerk.getY(),
1118 3.0e-12 * keplerianJerk.getNorm());
1119 Assertions.assertEquals(differentiate(orbit, shifted -> shifted.getPVCoordinates().getAcceleration().getZ()),
1120 keplerianJerk.getZ(),
1121 3.0e-12 * keplerianJerk.getNorm());
1122
1123 Assertions.assertTrue(Double.isNaN(orbit.getADot()));
1124 Assertions.assertTrue(Double.isNaN(orbit.getEquinoctialExDot()));
1125 Assertions.assertTrue(Double.isNaN(orbit.getEquinoctialEyDot()));
1126 Assertions.assertTrue(Double.isNaN(orbit.getHxDot()));
1127 Assertions.assertTrue(Double.isNaN(orbit.getHyDot()));
1128 Assertions.assertTrue(Double.isNaN(orbit.getLvDot()));
1129 Assertions.assertTrue(Double.isNaN(orbit.getLEDot()));
1130 Assertions.assertTrue(Double.isNaN(orbit.getLMDot()));
1131 Assertions.assertTrue(Double.isNaN(orbit.getEDot()));
1132 Assertions.assertTrue(Double.isNaN(orbit.getIDot()));
1133 Assertions.assertTrue(Double.isNaN(orbit.getPerigeeArgumentDot()));
1134 Assertions.assertTrue(Double.isNaN(orbit.getRightAscensionOfAscendingNodeDot()));
1135 Assertions.assertTrue(Double.isNaN(orbit.getTrueAnomalyDot()));
1136 Assertions.assertTrue(Double.isNaN(orbit.getEccentricAnomalyDot()));
1137 Assertions.assertTrue(Double.isNaN(orbit.getMeanAnomalyDot()));
1138 Assertions.assertTrue(Double.isNaN(orbit.getAnomalyDot(PositionAngle.TRUE)));
1139 Assertions.assertTrue(Double.isNaN(orbit.getAnomalyDot(PositionAngle.ECCENTRIC)));
1140 Assertions.assertTrue(Double.isNaN(orbit.getAnomalyDot(PositionAngle.MEAN)));
1141
1142 }
1143
1144 private <S extends Function<KeplerianOrbit, Double>>
1145 double differentiate(KeplerianOrbit orbit, S picker) {
1146 final DSFactory factory = new DSFactory(1, 1);
1147 FiniteDifferencesDifferentiator differentiator = new FiniteDifferencesDifferentiator(8, 0.1);
1148 UnivariateDifferentiableFunction diff = differentiator.differentiate(new UnivariateFunction() {
1149 public double value(double dt) {
1150 return picker.apply(orbit.shiftedBy(dt));
1151 }
1152 });
1153 return diff.value(factory.variable(0, 0.0)).getPartialDerivative(1);
1154 }
1155
1156 @Test
1157 public void testNonKeplerianEllipticDerivatives() {
1158 final AbsoluteDate date = new AbsoluteDate("2003-05-01T00:00:20.000", TimeScalesFactory.getUTC());
1159 final Vector3D position = new Vector3D(6896874.444705, 1956581.072644, -147476.245054);
1160 final Vector3D velocity = new Vector3D(166.816407662, -1106.783301861, -7372.745712770);
1161 final Vector3D acceleration = new Vector3D(-7.466182457944, -2.118153357345, 0.160004048437);
1162 final TimeStampedPVCoordinates pv = new TimeStampedPVCoordinates(date, position, velocity, acceleration);
1163 final Frame frame = FramesFactory.getEME2000();
1164 final double mu = Constants.EIGEN5C_EARTH_MU;
1165 final KeplerianOrbit orbit = new KeplerianOrbit(pv, frame, mu);
1166
1167 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getA()),
1168 orbit.getADot(),
1169 4.3e-8);
1170 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getEquinoctialEx()),
1171 orbit.getEquinoctialExDot(),
1172 2.1e-15);
1173 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getEquinoctialEy()),
1174 orbit.getEquinoctialEyDot(),
1175 5.3e-16);
1176 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getHx()),
1177 orbit.getHxDot(),
1178 1.6e-15);
1179 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getHy()),
1180 orbit.getHyDot(),
1181 7.3e-17);
1182 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getLv()),
1183 orbit.getLvDot(),
1184 1.1e-14);
1185 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getLE()),
1186 orbit.getLEDot(),
1187 7.2e-15);
1188 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getLM()),
1189 orbit.getLMDot(),
1190 4.7e-15);
1191 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getE()),
1192 orbit.getEDot(),
1193 6.9e-16);
1194 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getI()),
1195 orbit.getIDot(),
1196 5.8e-16);
1197 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getPerigeeArgument()),
1198 orbit.getPerigeeArgumentDot(),
1199 1.5e-12);
1200 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getRightAscensionOfAscendingNode()),
1201 orbit.getRightAscensionOfAscendingNodeDot(),
1202 1.5e-15);
1203 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getTrueAnomaly()),
1204 orbit.getTrueAnomalyDot(),
1205 1.5e-12);
1206 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getEccentricAnomaly()),
1207 orbit.getEccentricAnomalyDot(),
1208 1.5e-12);
1209 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getMeanAnomaly()),
1210 orbit.getMeanAnomalyDot(),
1211 1.5e-12);
1212 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getAnomaly(PositionAngle.TRUE)),
1213 orbit.getAnomalyDot(PositionAngle.TRUE),
1214 1.5e-12);
1215 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getAnomaly(PositionAngle.ECCENTRIC)),
1216 orbit.getAnomalyDot(PositionAngle.ECCENTRIC),
1217 1.5e-12);
1218 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getAnomaly(PositionAngle.MEAN)),
1219 orbit.getAnomalyDot(PositionAngle.MEAN),
1220 1.5e-12);
1221
1222 }
1223
1224 @Test
1225 public void testNonKeplerianHyperbolicDerivatives() {
1226 final AbsoluteDate date = new AbsoluteDate("2003-05-01T00:00:20.000", TimeScalesFactory.getUTC());
1227 final Vector3D position = new Vector3D(224267911.905821, 290251613.109399, 45534292.777492);
1228 final Vector3D velocity = new Vector3D(-1494.068165293, 1124.771027677, 526.915286134);
1229 final Vector3D acceleration = new Vector3D(-0.001295920501, -0.002233045187, -0.000349906292);
1230 final TimeStampedPVCoordinates pv = new TimeStampedPVCoordinates(date, position, velocity, acceleration);
1231 final Frame frame = FramesFactory.getEME2000();
1232 final double mu = Constants.EIGEN5C_EARTH_MU;
1233 final KeplerianOrbit orbit = new KeplerianOrbit(pv, frame, mu);
1234
1235 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getA()),
1236 orbit.getADot(),
1237 9.6e-8);
1238 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getEquinoctialEx()),
1239 orbit.getEquinoctialExDot(),
1240 2.8e-16);
1241 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getEquinoctialEy()),
1242 orbit.getEquinoctialEyDot(),
1243 3.6e-15);
1244 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getHx()),
1245 orbit.getHxDot(),
1246 1.4e-15);
1247 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getHy()),
1248 orbit.getHyDot(),
1249 9.4e-16);
1250 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getLv()),
1251 orbit.getLvDot(),
1252 5.6e-16);
1253 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getLE()),
1254 orbit.getLEDot(),
1255 9.0e-16);
1256 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getLM()),
1257 orbit.getLMDot(),
1258 1.8e-15);
1259 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getE()),
1260 orbit.getEDot(),
1261 1.8e-15);
1262 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getI()),
1263 orbit.getIDot(),
1264 3.6e-15);
1265 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getPerigeeArgument()),
1266 orbit.getPerigeeArgumentDot(),
1267 9.4e-16);
1268 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getRightAscensionOfAscendingNode()),
1269 orbit.getRightAscensionOfAscendingNodeDot(),
1270 1.1e-15);
1271 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getTrueAnomaly()),
1272 orbit.getTrueAnomalyDot(),
1273 1.4e-15);
1274 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getEccentricAnomaly()),
1275 orbit.getEccentricAnomalyDot(),
1276 9.2e-16);
1277 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getMeanAnomaly()),
1278 orbit.getMeanAnomalyDot(),
1279 1.4e-15);
1280 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getAnomaly(PositionAngle.TRUE)),
1281 orbit.getAnomalyDot(PositionAngle.TRUE),
1282 1.4e-15);
1283 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getAnomaly(PositionAngle.ECCENTRIC)),
1284 orbit.getAnomalyDot(PositionAngle.ECCENTRIC),
1285 9.2e-16);
1286 Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getAnomaly(PositionAngle.MEAN)),
1287 orbit.getAnomalyDot(PositionAngle.MEAN),
1288 1.4e-15);
1289
1290 }
1291
1292 private <S extends Function<KeplerianOrbit, Double>>
1293 double differentiate(TimeStampedPVCoordinates pv, Frame frame, double mu, S picker) {
1294 final DSFactory factory = new DSFactory(1, 1);
1295 FiniteDifferencesDifferentiator differentiator = new FiniteDifferencesDifferentiator(8, 0.1);
1296 UnivariateDifferentiableFunction diff = differentiator.differentiate(new UnivariateFunction() {
1297 public double value(double dt) {
1298 return picker.apply(new KeplerianOrbit(pv.shiftedBy(dt), frame, mu));
1299 }
1300 });
1301 return diff.value(factory.variable(0, 0.0)).getPartialDerivative(1);
1302 }
1303
1304 @Test
1305 public void testPositionAngleDerivatives() {
1306 final AbsoluteDate date = new AbsoluteDate("2003-05-01T00:00:20.000", TimeScalesFactory.getUTC());
1307 final Vector3D position = new Vector3D(6896874.444705, 1956581.072644, -147476.245054);
1308 final Vector3D velocity = new Vector3D(166.816407662, -1106.783301861, -7372.745712770);
1309 final Vector3D acceleration = new Vector3D(-7.466182457944, -2.118153357345, 0.160004048437);
1310 final TimeStampedPVCoordinates pv = new TimeStampedPVCoordinates(date, position, velocity, acceleration);
1311 final Frame frame = FramesFactory.getEME2000();
1312 final double mu = Constants.EIGEN5C_EARTH_MU;
1313 final KeplerianOrbit orbit = new KeplerianOrbit(pv, frame, mu);
1314
1315 for (PositionAngle type : PositionAngle.values()) {
1316 final KeplerianOrbit rebuilt = new KeplerianOrbit(orbit.getA(),
1317 orbit.getE(),
1318 orbit.getI(),
1319 orbit.getPerigeeArgument(),
1320 orbit.getRightAscensionOfAscendingNode(),
1321 orbit.getAnomaly(type),
1322 orbit.getADot(),
1323 orbit.getEDot(),
1324 orbit.getIDot(),
1325 orbit.getPerigeeArgumentDot(),
1326 orbit.getRightAscensionOfAscendingNodeDot(),
1327 orbit.getAnomalyDot(type),
1328 type, orbit.getFrame(), orbit.getDate(), orbit.getMu());
1329 MatcherAssert.assertThat(rebuilt.getA(), relativelyCloseTo(orbit.getA(), 1));
1330 MatcherAssert.assertThat(rebuilt.getE(), relativelyCloseTo(orbit.getE(), 1));
1331 MatcherAssert.assertThat(rebuilt.getI(), relativelyCloseTo(orbit.getI(), 1));
1332 MatcherAssert.assertThat(rebuilt.getPerigeeArgument(), relativelyCloseTo(orbit.getPerigeeArgument(), 1));
1333 MatcherAssert.assertThat(rebuilt.getRightAscensionOfAscendingNode(), relativelyCloseTo(orbit.getRightAscensionOfAscendingNode(), 1));
1334 MatcherAssert.assertThat(rebuilt.getADot(), relativelyCloseTo(orbit.getADot(), 1));
1335 MatcherAssert.assertThat(rebuilt.getEDot(), relativelyCloseTo(orbit.getEDot(), 1));
1336 MatcherAssert.assertThat(rebuilt.getIDot(), relativelyCloseTo(orbit.getIDot(), 1));
1337 MatcherAssert.assertThat(rebuilt.getPerigeeArgumentDot(), relativelyCloseTo(orbit.getPerigeeArgumentDot(), 1));
1338 MatcherAssert.assertThat(rebuilt.getRightAscensionOfAscendingNodeDot(), relativelyCloseTo(orbit.getRightAscensionOfAscendingNodeDot(), 1));
1339 for (PositionAngle type2 : PositionAngle.values()) {
1340 MatcherAssert.assertThat(rebuilt.getAnomaly(type2), relativelyCloseTo(orbit.getAnomaly(type2), 1));
1341 MatcherAssert.assertThat(rebuilt.getAnomalyDot(type2), relativelyCloseTo(orbit.getAnomalyDot(type2), 1));
1342 }
1343 }
1344
1345 }
1346
1347 @Test
1348 public void testPositionAngleHyperbolicDerivatives() {
1349 final AbsoluteDate date = new AbsoluteDate("2003-05-01T00:00:20.000", TimeScalesFactory.getUTC());
1350 final Vector3D position = new Vector3D(224267911.905821, 290251613.109399, 45534292.777492);
1351 final Vector3D velocity = new Vector3D(-1494.068165293, 1124.771027677, 526.915286134);
1352 final Vector3D acceleration = new Vector3D(-0.001295920501, -0.002233045187, -0.000349906292);
1353 final TimeStampedPVCoordinates pv = new TimeStampedPVCoordinates(date, position, velocity, acceleration);
1354 final Frame frame = FramesFactory.getEME2000();
1355 final double mu = Constants.EIGEN5C_EARTH_MU;
1356 final KeplerianOrbit orbit = new KeplerianOrbit(pv, frame, mu);
1357
1358 for (PositionAngle type : PositionAngle.values()) {
1359 final KeplerianOrbit rebuilt = new KeplerianOrbit(orbit.getA(),
1360 orbit.getE(),
1361 orbit.getI(),
1362 orbit.getPerigeeArgument(),
1363 orbit.getRightAscensionOfAscendingNode(),
1364 orbit.getAnomaly(type),
1365 orbit.getADot(),
1366 orbit.getEDot(),
1367 orbit.getIDot(),
1368 orbit.getPerigeeArgumentDot(),
1369 orbit.getRightAscensionOfAscendingNodeDot(),
1370 orbit.getAnomalyDot(type),
1371 type, orbit.getFrame(), orbit.getDate(), orbit.getMu());
1372 MatcherAssert.assertThat(rebuilt.getA(), relativelyCloseTo(orbit.getA(), 1));
1373 MatcherAssert.assertThat(rebuilt.getE(), relativelyCloseTo(orbit.getE(), 1));
1374 MatcherAssert.assertThat(rebuilt.getI(), relativelyCloseTo(orbit.getI(), 1));
1375 MatcherAssert.assertThat(rebuilt.getPerigeeArgument(), relativelyCloseTo(orbit.getPerigeeArgument(), 1));
1376 MatcherAssert.assertThat(rebuilt.getRightAscensionOfAscendingNode(), relativelyCloseTo(orbit.getRightAscensionOfAscendingNode(), 1));
1377 MatcherAssert.assertThat(rebuilt.getADot(), relativelyCloseTo(orbit.getADot(), 1));
1378 MatcherAssert.assertThat(rebuilt.getEDot(), relativelyCloseTo(orbit.getEDot(), 1));
1379 MatcherAssert.assertThat(rebuilt.getIDot(), relativelyCloseTo(orbit.getIDot(), 1));
1380 MatcherAssert.assertThat(rebuilt.getPerigeeArgumentDot(), relativelyCloseTo(orbit.getPerigeeArgumentDot(), 1));
1381 MatcherAssert.assertThat(rebuilt.getRightAscensionOfAscendingNodeDot(), relativelyCloseTo(orbit.getRightAscensionOfAscendingNodeDot(), 1));
1382 for (PositionAngle type2 : PositionAngle.values()) {
1383 MatcherAssert.assertThat(rebuilt.getAnomaly(type2), relativelyCloseTo(orbit.getAnomaly(type2), 3));
1384 MatcherAssert.assertThat(rebuilt.getAnomalyDot(type2), relativelyCloseTo(orbit.getAnomalyDot(type2), 4));
1385 }
1386 }
1387
1388 }
1389
1390 @Test
1391 public void testSerialization()
1392 throws IOException, ClassNotFoundException {
1393 Vector3D position = new Vector3D(-29536113.0, 30329259.0, -100125.0);
1394 Vector3D velocity = new Vector3D(-2194.0, -2141.0, -8.0);
1395 PVCoordinates pvCoordinates = new PVCoordinates(position, velocity);
1396 KeplerianOrbit orbit = new KeplerianOrbit(pvCoordinates, FramesFactory.getEME2000(), date, mu);
1397 Assertions.assertEquals(42255170.003, orbit.getA(), 1.0e-3);
1398
1399 ByteArrayOutputStream bos = new ByteArrayOutputStream();
1400 ObjectOutputStream oos = new ObjectOutputStream(bos);
1401 oos.writeObject(orbit);
1402
1403 Assertions.assertTrue(bos.size() > 280);
1404 Assertions.assertTrue(bos.size() < 330);
1405
1406 ByteArrayInputStream bis = new ByteArrayInputStream(bos.toByteArray());
1407 ObjectInputStream ois = new ObjectInputStream(bis);
1408 KeplerianOrbit deserialized = (KeplerianOrbit) ois.readObject();
1409 Assertions.assertEquals(orbit.getA(), deserialized.getA(), 1.0e-10);
1410 Assertions.assertEquals(orbit.getE(), deserialized.getE(), 1.0e-10);
1411 Assertions.assertEquals(orbit.getI(), deserialized.getI(), 1.0e-10);
1412 Assertions.assertEquals(orbit.getPerigeeArgument(), deserialized.getPerigeeArgument(), 1.0e-10);
1413 Assertions.assertEquals(orbit.getRightAscensionOfAscendingNode(), deserialized.getRightAscensionOfAscendingNode(), 1.0e-10);
1414 Assertions.assertEquals(orbit.getTrueAnomaly(), deserialized.getTrueAnomaly(), 1.0e-10);
1415 Assertions.assertTrue(Double.isNaN(orbit.getADot()) && Double.isNaN(deserialized.getADot()));
1416 Assertions.assertTrue(Double.isNaN(orbit.getEDot()) && Double.isNaN(deserialized.getEDot()));
1417 Assertions.assertTrue(Double.isNaN(orbit.getIDot()) && Double.isNaN(deserialized.getIDot()));
1418 Assertions.assertTrue(Double.isNaN(orbit.getPerigeeArgumentDot()) && Double.isNaN(deserialized.getPerigeeArgumentDot()));
1419 Assertions.assertTrue(Double.isNaN(orbit.getRightAscensionOfAscendingNodeDot()) && Double.isNaN(deserialized.getRightAscensionOfAscendingNodeDot()));
1420 Assertions.assertTrue(Double.isNaN(orbit.getTrueAnomalyDot()) && Double.isNaN(deserialized.getTrueAnomalyDot()));
1421 Assertions.assertEquals(orbit.getDate(), deserialized.getDate());
1422 Assertions.assertEquals(orbit.getMu(), deserialized.getMu(), 1.0e-10);
1423 Assertions.assertEquals(orbit.getFrame().getName(), deserialized.getFrame().getName());
1424
1425 }
1426
1427 @Test
1428 public void testSerializationWithDerivatives()
1429 throws IOException, ClassNotFoundException {
1430 Vector3D position = new Vector3D(-29536113.0, 30329259.0, -100125.0);
1431 Vector3D velocity = new Vector3D(-2194.0, -2141.0, -8.0);
1432 double r2 = position.getNormSq();
1433 double r = FastMath.sqrt(r2);
1434 Vector3D acceleration = new Vector3D(-mu / (r * r2), position,
1435 1, new Vector3D(-0.1, 0.2, 0.3));
1436 PVCoordinates pvCoordinates = new PVCoordinates(position, velocity, acceleration);
1437 KeplerianOrbit orbit = new KeplerianOrbit(pvCoordinates, FramesFactory.getEME2000(), date, mu);
1438 Assertions.assertEquals(42255170.003, orbit.getA(), 1.0e-3);
1439
1440 ByteArrayOutputStream bos = new ByteArrayOutputStream();
1441 ObjectOutputStream oos = new ObjectOutputStream(bos);
1442 oos.writeObject(orbit);
1443
1444 Assertions.assertTrue(bos.size() > 330);
1445 Assertions.assertTrue(bos.size() < 380);
1446
1447 ByteArrayInputStream bis = new ByteArrayInputStream(bos.toByteArray());
1448 ObjectInputStream ois = new ObjectInputStream(bis);
1449 KeplerianOrbit deserialized = (KeplerianOrbit) ois.readObject();
1450 Assertions.assertEquals(orbit.getA(), deserialized.getA(), 1.0e-10);
1451 Assertions.assertEquals(orbit.getE(), deserialized.getE(), 1.0e-10);
1452 Assertions.assertEquals(orbit.getI(), deserialized.getI(), 1.0e-10);
1453 Assertions.assertEquals(orbit.getPerigeeArgument(), deserialized.getPerigeeArgument(), 1.0e-10);
1454 Assertions.assertEquals(orbit.getRightAscensionOfAscendingNode(), deserialized.getRightAscensionOfAscendingNode(), 1.0e-10);
1455 Assertions.assertEquals(orbit.getTrueAnomaly(), deserialized.getTrueAnomaly(), 1.0e-10);
1456 Assertions.assertEquals(orbit.getADot(), deserialized.getADot(), 1.0e-10);
1457 Assertions.assertEquals(orbit.getEDot(), deserialized.getEDot(), 1.0e-10);
1458 Assertions.assertEquals(orbit.getIDot(), deserialized.getIDot(), 1.0e-10);
1459 Assertions.assertEquals(orbit.getPerigeeArgumentDot(), deserialized.getPerigeeArgumentDot(), 1.0e-10);
1460 Assertions.assertEquals(orbit.getRightAscensionOfAscendingNodeDot(), deserialized.getRightAscensionOfAscendingNodeDot(), 1.0e-10);
1461 Assertions.assertEquals(orbit.getTrueAnomalyDot(), deserialized.getTrueAnomalyDot(), 1.0e-10);
1462 Assertions.assertEquals(orbit.getDate(), deserialized.getDate());
1463 Assertions.assertEquals(orbit.getMu(), deserialized.getMu(), 1.0e-10);
1464 Assertions.assertEquals(orbit.getFrame().getName(), deserialized.getFrame().getName());
1465
1466 }
1467
1468 @Test
1469 public void testEquatorialRetrograde() {
1470 Vector3D position = new Vector3D(10000000.0, 0.0, 0.0);
1471 Vector3D velocity = new Vector3D(0.0, -6500.0, 1.0e-10);
1472 double r2 = position.getNormSq();
1473 double r = FastMath.sqrt(r2);
1474 Vector3D acceleration = new Vector3D(-mu / (r * r2), position,
1475 1, new Vector3D(-0.1, 0.2, 0.3));
1476 PVCoordinates pvCoordinates = new PVCoordinates(position, velocity, acceleration);
1477 KeplerianOrbit orbit = new KeplerianOrbit(pvCoordinates, FramesFactory.getEME2000(), date, mu);
1478 Assertions.assertEquals(10637829.465, orbit.getA(), 1.0e-3);
1479 Assertions.assertEquals(-738.145, orbit.getADot(), 1.0e-3);
1480 Assertions.assertEquals(0.05995861, orbit.getE(), 1.0e-8);
1481 Assertions.assertEquals(-6.523e-5, orbit.getEDot(), 1.0e-8);
1482 Assertions.assertEquals(FastMath.PI, orbit.getI(), 2.0e-14);
1483 Assertions.assertEquals(-4.615e-5, orbit.getIDot(), 1.0e-8);
1484 Assertions.assertTrue(Double.isNaN(orbit.getHx()));
1485 Assertions.assertTrue(Double.isNaN(orbit.getHxDot()));
1486 Assertions.assertTrue(Double.isNaN(orbit.getHy()));
1487 Assertions.assertTrue(Double.isNaN(orbit.getHyDot()));
1488 }
1489
1490 @Test
1491 public void testDerivativesConversionSymmetry() {
1492 final AbsoluteDate date = new AbsoluteDate("2003-05-01T00:01:20.000", TimeScalesFactory.getUTC());
1493 Vector3D position = new Vector3D(6893443.400234382, 1886406.1073757345, -589265.1150359757);
1494 Vector3D velocity = new Vector3D(-281.1261461082365, -1231.6165642450928, -7348.756363469432);
1495 Vector3D acceleration = new Vector3D(-7.460341170581685, -2.0415957334584527, 0.6393322823627762);
1496 PVCoordinates pvCoordinates = new PVCoordinates( position, velocity, acceleration);
1497 KeplerianOrbit orbit = new KeplerianOrbit(pvCoordinates, FramesFactory.getEME2000(),
1498 date, Constants.EIGEN5C_EARTH_MU);
1499 Assertions.assertTrue(orbit.hasDerivatives());
1500 double r2 = position.getNormSq();
1501 double r = FastMath.sqrt(r2);
1502 Vector3D keplerianAcceleration = new Vector3D(-orbit.getMu() / (r2 * r), position);
1503 Assertions.assertEquals(0.0101, Vector3D.distance(keplerianAcceleration, acceleration), 1.0e-4);
1504
1505 for (OrbitType type : OrbitType.values()) {
1506 Orbit converted = type.convertType(orbit);
1507 Assertions.assertTrue(converted.hasDerivatives());
1508 KeplerianOrbit rebuilt = (KeplerianOrbit) OrbitType.KEPLERIAN.convertType(converted);
1509 Assertions.assertTrue(rebuilt.hasDerivatives());
1510 Assertions.assertEquals(orbit.getADot(), rebuilt.getADot(), 3.0e-13);
1511 Assertions.assertEquals(orbit.getEDot(), rebuilt.getEDot(), 1.0e-15);
1512 Assertions.assertEquals(orbit.getIDot(), rebuilt.getIDot(), 1.0e-15);
1513 Assertions.assertEquals(orbit.getPerigeeArgumentDot(), rebuilt.getPerigeeArgumentDot(), 1.0e-15);
1514 Assertions.assertEquals(orbit.getRightAscensionOfAscendingNodeDot(), rebuilt.getRightAscensionOfAscendingNodeDot(), 1.0e-15);
1515 Assertions.assertEquals(orbit.getTrueAnomalyDot(), rebuilt.getTrueAnomalyDot(), 1.0e-15);
1516 }
1517
1518 }
1519
1520 @Test
1521 public void testDerivativesConversionSymmetryHyperbolic() {
1522 final AbsoluteDate date = new AbsoluteDate("2003-05-01T00:00:20.000", TimeScalesFactory.getUTC());
1523 final Vector3D position = new Vector3D(224267911.905821, 290251613.109399, 45534292.777492);
1524 final Vector3D velocity = new Vector3D(-1494.068165293, 1124.771027677, 526.915286134);
1525 final Vector3D acceleration = new Vector3D(-0.001295920501, -0.002233045187, -0.000349906292);
1526 PVCoordinates pvCoordinates = new PVCoordinates( position, velocity, acceleration);
1527 KeplerianOrbit orbit = new KeplerianOrbit(pvCoordinates, FramesFactory.getEME2000(),
1528 date, Constants.EIGEN5C_EARTH_MU);
1529 Assertions.assertTrue(orbit.hasDerivatives());
1530 double r2 = position.getNormSq();
1531 double r = FastMath.sqrt(r2);
1532 Vector3D keplerianAcceleration = new Vector3D(-orbit.getMu() / (r2 * r), position);
1533 Assertions.assertEquals(4.78e-4, Vector3D.distance(keplerianAcceleration, acceleration), 1.0e-6);
1534
1535 OrbitType type = OrbitType.CARTESIAN;
1536 Orbit converted = type.convertType(orbit);
1537 Assertions.assertTrue(converted.hasDerivatives());
1538 KeplerianOrbit rebuilt = (KeplerianOrbit) OrbitType.KEPLERIAN.convertType(converted);
1539 Assertions.assertTrue(rebuilt.hasDerivatives());
1540 Assertions.assertEquals(orbit.getADot(), rebuilt.getADot(), 3.0e-13);
1541 Assertions.assertEquals(orbit.getEDot(), rebuilt.getEDot(), 1.0e-15);
1542 Assertions.assertEquals(orbit.getIDot(), rebuilt.getIDot(), 1.0e-15);
1543 Assertions.assertEquals(orbit.getPerigeeArgumentDot(), rebuilt.getPerigeeArgumentDot(), 1.0e-15);
1544 Assertions.assertEquals(orbit.getRightAscensionOfAscendingNodeDot(), rebuilt.getRightAscensionOfAscendingNodeDot(), 1.0e-15);
1545 Assertions.assertEquals(orbit.getTrueAnomalyDot(), rebuilt.getTrueAnomalyDot(), 1.0e-15);
1546
1547 }
1548
1549 @Test
1550 public void testToString() {
1551 Vector3D position = new Vector3D(-29536113.0, 30329259.0, -100125.0);
1552 Vector3D velocity = new Vector3D(-2194.0, -2141.0, -8.0);
1553 PVCoordinates pvCoordinates = new PVCoordinates(position, velocity);
1554 KeplerianOrbit orbit = new KeplerianOrbit(pvCoordinates, FramesFactory.getEME2000(), date, mu);
1555 Assertions.assertEquals("Keplerian parameters: {a: 4.225517000282565E7; e: 0.002146216321416967; i: 0.20189257051515358; pa: 13.949966363606599; raan: -87.91788415673473; v: -151.79096272977213;}",
1556 orbit.toString());
1557 }
1558
1559 @Test
1560 public void testCopyNonKeplerianAcceleration() {
1561
1562 final Frame eme2000 = FramesFactory.getEME2000();
1563
1564
1565 final Vector3D position = new Vector3D(42164140, 0, 0);
1566
1567 final PVCoordinates pv = new PVCoordinates(position,
1568 new Vector3D(0, FastMath.sqrt(mu / position.getNorm()), 0));
1569
1570 final Orbit orbit = new KeplerianOrbit(pv, eme2000, date, mu);
1571
1572
1573 final Orbit orbitCopy = new KeplerianOrbit(orbit);
1574
1575
1576 final Orbit shiftedOrbit = orbit.shiftedBy(10);
1577 final Orbit shiftedOrbitCopy = orbitCopy.shiftedBy(10);
1578
1579 Assertions.assertEquals(0.0,
1580 Vector3D.distance(shiftedOrbit.getPVCoordinates().getPosition(),
1581 shiftedOrbitCopy.getPVCoordinates().getPosition()),
1582 1.0e-10);
1583 Assertions.assertEquals(0.0,
1584 Vector3D.distance(shiftedOrbit.getPVCoordinates().getVelocity(),
1585 shiftedOrbitCopy.getPVCoordinates().getVelocity()),
1586 1.0e-10);
1587
1588 }
1589
1590 @Test
1591 public void testIssue674() {
1592 try {
1593 new KeplerianOrbit(24464560.0, -0.7311, 0.122138, 3.10686, 1.00681,
1594 0.048363, PositionAngle.MEAN,
1595 FramesFactory.getEME2000(), date, mu);
1596 Assertions.fail("an exception should have been thrown");
1597 } catch (OrekitException oe) {
1598 Assertions.assertEquals(OrekitMessages.INVALID_PARAMETER_RANGE, oe.getSpecifier());
1599 Assertions.assertEquals("eccentricity", oe.getParts()[0]);
1600 Assertions.assertEquals(-0.7311, ((Double) oe.getParts()[1]).doubleValue(), Double.MIN_VALUE);
1601 Assertions.assertEquals(0.0, ((Double) oe.getParts()[2]).doubleValue(), Double.MIN_VALUE);
1602 Assertions.assertEquals(Double.POSITIVE_INFINITY, ((Double) oe.getParts()[3]).doubleValue(), Double.MIN_VALUE);
1603 }
1604 }
1605
1606 @Test
1607 public void testNormalize() {
1608 KeplerianOrbit withoutDerivatives =
1609 new KeplerianOrbit(42166712.0, 0.005, 1.6,
1610 -0.3, 1.25, 0.4, PositionAngle.MEAN,
1611 FramesFactory.getEME2000(), date, mu);
1612 KeplerianOrbit ref =
1613 new KeplerianOrbit(24000000.0, 0.012, 1.9,
1614 -6.28, 6.28, 12.56, PositionAngle.MEAN,
1615 FramesFactory.getEME2000(), date, mu);
1616
1617 KeplerianOrbit normalized1 = (KeplerianOrbit) OrbitType.KEPLERIAN.normalize(withoutDerivatives, ref);
1618 Assertions.assertFalse(normalized1.hasDerivatives());
1619 Assertions.assertEquals(0.0, normalized1.getA() - withoutDerivatives.getA(), 1.0e-6);
1620 Assertions.assertEquals(0.0, normalized1.getE() - withoutDerivatives.getE(), 1.0e-10);
1621 Assertions.assertEquals(0.0, normalized1.getI() - withoutDerivatives.getI(), 1.0e-10);
1622 Assertions.assertEquals(-MathUtils.TWO_PI, normalized1.getPerigeeArgument() - withoutDerivatives.getPerigeeArgument(), 1.0e-10);
1623 Assertions.assertEquals(+MathUtils.TWO_PI, normalized1.getRightAscensionOfAscendingNode() - withoutDerivatives.getRightAscensionOfAscendingNode(), 1.0e-10);
1624 Assertions.assertEquals(2 * MathUtils.TWO_PI, normalized1.getTrueAnomaly() - withoutDerivatives.getTrueAnomaly(), 1.0e-10);
1625 Assertions.assertTrue(Double.isNaN(normalized1.getADot()));
1626 Assertions.assertTrue(Double.isNaN(normalized1.getEDot()));
1627 Assertions.assertTrue(Double.isNaN(normalized1.getIDot()));
1628 Assertions.assertTrue(Double.isNaN(normalized1.getPerigeeArgumentDot()));
1629 Assertions.assertTrue(Double.isNaN(normalized1.getRightAscensionOfAscendingNodeDot()));
1630 Assertions.assertTrue(Double.isNaN(normalized1.getTrueAnomalyDot()));
1631
1632 double[] p = new double[6];
1633 OrbitType.KEPLERIAN.mapOrbitToArray(withoutDerivatives, PositionAngle.TRUE, p, null);
1634 KeplerianOrbit withDerivatives = (KeplerianOrbit) OrbitType.KEPLERIAN.mapArrayToOrbit(p,
1635 new double[] { 1.0, 2.0, 3.0, 4.0, 5.0, 6.0 },
1636 PositionAngle.TRUE,
1637 withoutDerivatives.getDate(),
1638 withoutDerivatives.getMu(),
1639 withoutDerivatives.getFrame());
1640 KeplerianOrbit normalized2 = (KeplerianOrbit) OrbitType.KEPLERIAN.normalize(withDerivatives, ref);
1641 Assertions.assertTrue(normalized2.hasDerivatives());
1642 Assertions.assertEquals(0.0, normalized2.getA() - withoutDerivatives.getA(), 1.0e-6);
1643 Assertions.assertEquals(0.0, normalized2.getE() - withoutDerivatives.getE(), 1.0e-10);
1644 Assertions.assertEquals(0.0, normalized2.getI() - withoutDerivatives.getI(), 1.0e-10);
1645 Assertions.assertEquals(-MathUtils.TWO_PI, normalized2.getPerigeeArgument() - withoutDerivatives.getPerigeeArgument(), 1.0e-10);
1646 Assertions.assertEquals(+MathUtils.TWO_PI, normalized2.getRightAscensionOfAscendingNode() - withoutDerivatives.getRightAscensionOfAscendingNode(), 1.0e-10);
1647 Assertions.assertEquals(2 * MathUtils.TWO_PI, normalized2.getTrueAnomaly() - withoutDerivatives.getTrueAnomaly(), 1.0e-10);
1648 Assertions.assertEquals(0.0, normalized2.getADot() - withDerivatives.getADot(), 1.0e-10);
1649 Assertions.assertEquals(0.0, normalized2.getEDot() - withDerivatives.getEDot(), 1.0e-10);
1650 Assertions.assertEquals(0.0, normalized2.getIDot() - withDerivatives.getIDot(), 1.0e-10);
1651 Assertions.assertEquals(0.0, normalized2.getPerigeeArgumentDot() - withDerivatives.getPerigeeArgumentDot(), 1.0e-10);
1652 Assertions.assertEquals(0.0, normalized2.getRightAscensionOfAscendingNodeDot() - withDerivatives.getRightAscensionOfAscendingNodeDot(), 1.0e-10);
1653 Assertions.assertEquals(0.0, normalized2.getTrueAnomalyDot() - withDerivatives.getTrueAnomalyDot(), 1.0e-10);
1654
1655 }
1656
1657 @Test
1658 public void testKeplerianToPvToKeplerian() {
1659
1660 Frame eci = FramesFactory.getGCRF();
1661 AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
1662 double mu = Constants.EGM96_EARTH_MU;
1663 KeplerianOrbit orbit = new KeplerianOrbit(6878e3, 0, 0, 0, 0, 0,
1664 PositionAngle.TRUE, eci, date, mu);
1665
1666
1667 KeplerianOrbit pvOrbit = new KeplerianOrbit(orbit.getPVCoordinates(), eci, mu);
1668
1669
1670 Assertions.assertEquals(orbit.shiftedBy(1).getA(), pvOrbit.shiftedBy(1).getA(), 1e-6);
1671 Assertions.assertEquals(orbit.shiftedBy(1).getE(), pvOrbit.shiftedBy(1).getE(), 1e-14);
1672 Assertions.assertEquals(orbit.shiftedBy(1).getI(), pvOrbit.shiftedBy(1).getI(), 1e-10);
1673 Assertions.assertEquals(orbit.shiftedBy(1).getPerigeeArgument(), pvOrbit.shiftedBy(1).getPerigeeArgument(), 1e-10);
1674 Assertions.assertEquals(orbit.shiftedBy(1).getRightAscensionOfAscendingNode(), pvOrbit.shiftedBy(1).getRightAscensionOfAscendingNode(), 1e-10);
1675 Assertions.assertEquals(orbit.shiftedBy(1).getTrueAnomaly(), pvOrbit.shiftedBy(1).getTrueAnomaly(), 1e-10);
1676
1677 }
1678
1679
1680 @BeforeEach
1681 public void setUp() {
1682
1683 Utils.setDataRoot("regular-data");
1684
1685
1686 date = AbsoluteDate.J2000_EPOCH;
1687
1688
1689 mu = 3.9860047e14;
1690 }
1691
1692 @AfterEach
1693 public void tearDown() {
1694 date = null;
1695 }
1696
1697 }