1   /* Copyright 2002-2025 CS GROUP
2    * Licensed to CS GROUP (CS) under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * CS licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *   http://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  package org.orekit.frames;
18  
19  import java.io.IOException;
20  import java.util.stream.Stream;
21  
22  import org.hipparchus.CalculusFieldElement;
23  import org.hipparchus.Field;
24  import org.hipparchus.analysis.polynomials.PolynomialFunction;
25  import org.hipparchus.complex.Complex;
26  import org.hipparchus.complex.ComplexField;
27  import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
28  import org.hipparchus.geometry.euclidean.threed.Vector3D;
29  import org.hipparchus.random.RandomGenerator;
30  import org.hipparchus.random.Well1024a;
31  import org.hipparchus.util.Binary64Field;
32  import org.hipparchus.util.Binary64;
33  import org.hipparchus.util.FastMath;
34  import org.hipparchus.util.MathUtils;
35  import org.hipparchus.util.SinCos;
36  import org.hipparchus.util.FieldSinCos;
37  import org.junit.jupiter.api.AfterEach;
38  import org.junit.jupiter.api.Assertions;
39  import org.junit.jupiter.api.BeforeEach;
40  import org.junit.jupiter.api.Test;
41  import org.junit.jupiter.params.ParameterizedTest;
42  import org.junit.jupiter.params.provider.Arguments;
43  import org.junit.jupiter.params.provider.MethodSource;
44  import org.orekit.Utils;
45  import org.orekit.bodies.BodyShape;
46  import org.orekit.bodies.FieldGeodeticPoint;
47  import org.orekit.bodies.GeodeticPoint;
48  import org.orekit.bodies.OneAxisEllipsoid;
49  import org.orekit.orbits.CircularOrbit;
50  import org.orekit.orbits.FieldCircularOrbit;
51  import org.orekit.orbits.PositionAngleType;
52  import org.orekit.propagation.FieldSpacecraftState;
53  import org.orekit.propagation.SpacecraftState;
54  import org.orekit.propagation.analytical.FieldKeplerianPropagator;
55  import org.orekit.propagation.analytical.KeplerianPropagator;
56  import org.orekit.time.AbsoluteDate;
57  import org.orekit.time.DateComponents;
58  import org.orekit.time.FieldAbsoluteDate;
59  import org.orekit.time.TimeComponents;
60  import org.orekit.time.TimeScalesFactory;
61  import org.orekit.utils.Constants;
62  import org.orekit.utils.FieldPVCoordinates;
63  import org.orekit.utils.FieldTrackingCoordinates;
64  import org.orekit.utils.IERSConventions;
65  import org.orekit.utils.PVCoordinates;
66  import org.orekit.utils.TrackingCoordinates;
67  
68  
69  class TopocentricFrameTest {
70  
71      // Computation date
72      private AbsoluteDate date;
73  
74      // Reference frame = ITRF
75      private Frame itrf;
76  
77      // Earth shape
78      OneAxisEllipsoid earthSpheric;
79  
80      // Body mu
81      private double mu;
82  
83  
84      @Test
85      void testZero() {
86  
87          final GeodeticPoint point = new GeodeticPoint(0., 0., 0.);
88          final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, point, "zero");
89  
90          // Check that frame directions are aligned
91          final double xDiff = Vector3D.dotProduct(topoFrame.getEast(), Vector3D.PLUS_J);
92          final double yDiff = Vector3D.dotProduct(topoFrame.getNorth(), Vector3D.PLUS_K);
93          final double zDiff = Vector3D.dotProduct(topoFrame.getZenith(), Vector3D.PLUS_I);
94          Assertions.assertEquals(1., xDiff, Utils.epsilonTest);
95          Assertions.assertEquals(1., yDiff, Utils.epsilonTest);
96          Assertions.assertEquals(1., zDiff, Utils.epsilonTest);
97     }
98  
99      @Test
100     void testPole() {
101 
102         final GeodeticPoint point = new GeodeticPoint(FastMath.PI/2., 0., 0.);
103         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, point, "north pole");
104 
105         // Check that frame directions are aligned
106         final double xDiff = Vector3D.dotProduct(topoFrame.getEast(), Vector3D.PLUS_J);
107         final double yDiff = Vector3D.dotProduct(topoFrame.getSouth(), Vector3D.PLUS_I);
108         final double zDiff = Vector3D.dotProduct(topoFrame.getZenith(), Vector3D.PLUS_K);
109         Assertions.assertEquals(1., xDiff, Utils.epsilonTest);
110         Assertions.assertEquals(1., yDiff, Utils.epsilonTest);
111         Assertions.assertEquals(1., zDiff, Utils.epsilonTest);
112    }
113 
114     @Test
115     void testNormalLatitudes() {
116 
117         // First point at latitude 45°
118         final GeodeticPoint point1 = new GeodeticPoint(FastMath.toRadians(45.), FastMath.toRadians(30.), 0.);
119         final TopocentricFrame topoFrame1 = new TopocentricFrame(earthSpheric, point1, "lat 45");
120         Assertions.assertEquals(0.0,
121                                 Vector3D.distance(topoFrame1.getCartesianPoint(), earthSpheric.transform(point1)),
122                                 7.0e-9);
123 
124         // Second point at latitude -45° and same longitude
125         final GeodeticPoint point2 = new GeodeticPoint(FastMath.toRadians(-45.), FastMath.toRadians(30.), 0.);
126         final TopocentricFrame topoFrame2 = new TopocentricFrame(earthSpheric, point2, "lat -45");
127         Assertions.assertEquals(0.0,
128                                 Vector3D.distance(topoFrame2.getCartesianPoint(), earthSpheric.transform(point2)),
129                                 2.0e-9);
130 
131         // Check that frame North and Zenith directions are all normal to each other, and East are the same
132         final double xDiff = Vector3D.dotProduct(topoFrame1.getEast(), topoFrame2.getEast());
133         final double yDiff = Vector3D.dotProduct(topoFrame1.getNorth(), topoFrame2.getNorth());
134         final double zDiff = Vector3D.dotProduct(topoFrame1.getZenith(), topoFrame2.getZenith());
135 
136         Assertions.assertEquals(1., xDiff, Utils.epsilonTest);
137         Assertions.assertEquals(0., yDiff, Utils.epsilonTest);
138         Assertions.assertEquals(0., zDiff, Utils.epsilonTest);
139   }
140 
141     @Test
142     void testOppositeLongitudes() {
143 
144         // First point at latitude 45°
145         final GeodeticPoint point1 = new GeodeticPoint(FastMath.toRadians(45.), FastMath.toRadians(30.), 0.);
146         final TopocentricFrame topoFrame1 = new TopocentricFrame(earthSpheric, point1, "lon 30");
147         final GeodeticPoint p1 = topoFrame1.getPoint();
148         Assertions.assertEquals(point1.getLatitude(), p1.getLatitude(), 1.0e-15);
149         Assertions.assertEquals(point1.getLongitude(), p1.getLongitude(), 1.0e-15);
150         Assertions.assertEquals(point1.getAltitude(), p1.getAltitude(), 1.0e-15);
151 
152         // Second point at latitude -45° and same longitude
153         final GeodeticPoint point2 = new GeodeticPoint(FastMath.toRadians(45.), FastMath.toRadians(210.), 0.);
154         final TopocentricFrame topoFrame2 = new TopocentricFrame(earthSpheric, point2, "lon 210");
155 
156         // Check that frame North and Zenith directions are all normal to each other,
157         // and East of the one is West of the other
158         final double xDiff = Vector3D.dotProduct(topoFrame1.getEast(), topoFrame2.getWest());
159         final double yDiff = Vector3D.dotProduct(topoFrame1.getNorth(), topoFrame2.getNorth());
160         final double zDiff = Vector3D.dotProduct(topoFrame1.getZenith(), topoFrame2.getZenith());
161 
162         Assertions.assertEquals(1., xDiff, Utils.epsilonTest);
163         Assertions.assertEquals(0., yDiff, Utils.epsilonTest);
164         Assertions.assertEquals(0., zDiff, Utils.epsilonTest);
165   }
166 
167     @Test
168     void testAntipodes() {
169 
170         // First point at latitude 45° and longitude 30
171         final GeodeticPoint point1 = new GeodeticPoint(FastMath.toRadians(45.), FastMath.toRadians(30.), 0.);
172         final TopocentricFrame topoFrame1 = new TopocentricFrame(earthSpheric, point1, "lon 30");
173 
174         // Second point at latitude -45° and longitude 210
175         final GeodeticPoint point2 = new GeodeticPoint(FastMath.toRadians(-45.), FastMath.toRadians(210.), 0.);
176         final TopocentricFrame topoFrame2 = new TopocentricFrame(earthSpheric, point2, "lon 210");
177 
178         // Check that frame Zenith directions are opposite to each other,
179         // and East and North are the same
180         final double xDiff = Vector3D.dotProduct(topoFrame1.getEast(), topoFrame2.getWest());
181         final double yDiff = Vector3D.dotProduct(topoFrame1.getNorth(), topoFrame2.getNorth());
182         final double zDiff = Vector3D.dotProduct(topoFrame1.getZenith(), topoFrame2.getZenith());
183 
184         Assertions.assertEquals(1., xDiff, Utils.epsilonTest);
185         Assertions.assertEquals(1., yDiff, Utils.epsilonTest);
186         Assertions.assertEquals(-1., zDiff, Utils.epsilonTest);
187 
188         Assertions.assertEquals(1, Vector3D.dotProduct(topoFrame1.getNadir(), topoFrame2.getZenith()), Utils.epsilonTest);
189         Assertions.assertEquals(1, Vector3D.dotProduct(topoFrame1.getZenith(), topoFrame2.getNadir()), Utils.epsilonTest);
190 
191     }
192 
193     @Test
194     void testSiteAtZenith() {
195 
196         // Surface point at latitude 45°
197         final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(45.), FastMath.toRadians(30.), 0.);
198         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, point, "lon 30 lat 45");
199 
200         // Point at 800 km over zenith
201         final GeodeticPoint satPoint = new GeodeticPoint(FastMath.toRadians(45.), FastMath.toRadians(30.), 800000.);
202         final TrackingCoordinates tc = topoFrame.getTrackingCoordinates(earthSpheric.transform(satPoint),
203                                                                         earthSpheric.getBodyFrame(),
204                                                                         date);
205 
206         // Zenith point elevation = 90 deg
207         Assertions.assertEquals(FastMath.PI/2., tc.getElevation(), Utils.epsilonAngle);
208 
209         // Zenith point range = defined altitude
210         Assertions.assertEquals(800000., tc.getRange(), 1e-8);
211   }
212 
213     @Test
214     void testFieldSiteAtZenith() {
215         doTestFieldSiteAtZenith(Binary64Field.getInstance());
216     }
217 
218     private <T extends CalculusFieldElement<T>> void doTestFieldSiteAtZenith(final Field<T> field) {
219 
220         // zero
221         final T zero = field.getZero();
222 
223         // Surface point at latitude 45°
224         final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(45.), FastMath.toRadians(30.), 0.);
225         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, point, "lon 30 lat 45");
226 
227         // Point at 800 km over zenith
228         final FieldGeodeticPoint<T> satPoint = new FieldGeodeticPoint<>(zero.add(FastMath.toRadians(45.)),
229                                                                         zero.add(FastMath.toRadians(30.)),
230                                                                         zero.add(800000.));
231         final FieldAbsoluteDate<T> fieldDate = new FieldAbsoluteDate<>(field, date);
232         final FieldTrackingCoordinates<T> tc = topoFrame.getTrackingCoordinates(earthSpheric.transform(satPoint),
233                                                                                 earthSpheric.getBodyFrame(),
234                                                                                 fieldDate);
235 
236         // Zenith point elevation = 90 deg
237         Assertions.assertEquals(FastMath.PI/2., tc.getElevation().getReal(), Utils.epsilonAngle);
238 
239         // Zenith point range = defined altitude
240         Assertions.assertEquals(800000., tc.getRange().getReal(), 1e-8);
241 
242     }
243 
244     @Test
245     void testAzimuthEquatorial() {
246 
247         // Surface point at latitude 0
248         final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(0.), FastMath.toRadians(30.), 0.);
249         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, point, "lon 30 lat 0");
250 
251         // Point at infinite, separated by +20 deg in longitude
252         // *****************************************************
253         GeodeticPoint infPoint = new GeodeticPoint(FastMath.toRadians(0.), FastMath.toRadians(50.), 1000000000.);
254         TrackingCoordinates tc = topoFrame.getTrackingCoordinates(earthSpheric.transform(infPoint), earthSpheric.getBodyFrame(), date);
255 
256         // Azimuth = pi/2
257         Assertions.assertEquals(FastMath.PI/2., tc.getAzimuth(), Utils.epsilonAngle);
258 
259         // Elevation = pi/2 - longitude difference
260         Assertions.assertEquals(FastMath.PI/2. - FastMath.abs(point.getLongitude() - infPoint.getLongitude()),
261                                 tc.getElevation(),
262                                 1.e-2);
263 
264         // Point at infinite, separated by -20 deg in longitude
265         // *****************************************************
266         infPoint = new GeodeticPoint(FastMath.toRadians(0.), FastMath.toRadians(10.), 1000000000.);
267         tc = topoFrame.getTrackingCoordinates(earthSpheric.transform(infPoint), earthSpheric.getBodyFrame(), date);
268 
269         // Azimuth = pi/2
270         Assertions.assertEquals(3*FastMath.PI/2., tc.getAzimuth(), Utils.epsilonAngle);
271 
272         // Elevation = pi/2 - longitude difference
273         Assertions.assertEquals(FastMath.PI/2. - FastMath.abs(point.getLongitude() - infPoint.getLongitude()),
274                                 tc.getElevation(),
275                                 1.e-2);
276 
277     }
278 
279     @Test
280     void testFieldAzimuthEquatorial() {
281         doTestFieldAzimuthEquatorial(Binary64Field.getInstance());
282     }
283 
284     private <T extends CalculusFieldElement<T>> void doTestFieldAzimuthEquatorial(final Field<T> field) {
285 
286         // zero
287         final T zero = field.getZero();
288 
289         // Surface point at latitude 0
290         final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(0.), FastMath.toRadians(30.), 0.);
291         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, point, "lon 30 lat 0");
292 
293         // Point at infinite, separated by +20 deg in longitude
294         // *****************************************************
295         FieldGeodeticPoint<T> infPoint = new FieldGeodeticPoint<>(zero.add(FastMath.toRadians(0.)),
296                                                                   zero.add(FastMath.toRadians(50.)),
297                                                                   zero.add(1000000000.));
298 
299         // Field date
300         final FieldAbsoluteDate<T> fieldDate = new FieldAbsoluteDate<>(field, date);
301         FieldTrackingCoordinates<T> tc = topoFrame.getTrackingCoordinates(earthSpheric.transform(infPoint),
302                                                                           earthSpheric.getBodyFrame(),
303                                                                           fieldDate);
304         // Azimuth = pi/2
305         Assertions.assertEquals(FastMath.PI/2., tc.getAzimuth().getReal(), Utils.epsilonAngle);
306 
307         // Elevation = pi/2 - longitude difference
308         Assertions.assertEquals(FastMath.abs(infPoint.getLongitude().negate().add(point.getLongitude())).negate().add(FastMath.PI/2.).getReal(),
309                                 tc.getElevation().getReal(),
310                                 1.e-2);
311 
312         // Point at infinite, separated by -20 deg in longitude
313         // *****************************************************
314         infPoint = new FieldGeodeticPoint<>(zero.add(FastMath.toRadians(0.)),
315                                             zero.add(FastMath.toRadians(10.)),
316                                             zero.add(1000000000.));
317 
318         // Azimuth = pi/2
319         tc = topoFrame.getTrackingCoordinates(earthSpheric.transform(infPoint), earthSpheric.getBodyFrame(), fieldDate);
320         Assertions.assertEquals(3*FastMath.PI/2., tc.getAzimuth().getReal(), Utils.epsilonAngle);
321 
322         // Site = pi/2 - longitude difference
323         Assertions.assertEquals(FastMath.abs(infPoint.getLongitude().negate().add(point.getLongitude())).negate().add(FastMath.PI/2.).getReal(),
324                                 tc.getElevation().getReal(),
325                                 1.e-2);
326 
327     }
328 
329     @Test
330     void testAzimuthPole() {
331 
332         // Surface point at latitude 0
333         final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(89.999), FastMath.toRadians(0.), 0.);
334         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, point, "lon 0 lat 90");
335 
336         // Point at 30 deg longitude
337         // **************************
338         GeodeticPoint satPoint = new GeodeticPoint(FastMath.toRadians(28.), FastMath.toRadians(30.), 800000.);
339         TrackingCoordinates tc = topoFrame.getTrackingCoordinates(earthSpheric.transform(satPoint), earthSpheric.getBodyFrame(), date);
340 
341         Assertions.assertEquals(FastMath.PI - satPoint.getLongitude(), tc.getAzimuth(), 1.e-5);
342 
343         // Point at -30 deg longitude
344         // ***************************
345         satPoint = new GeodeticPoint(FastMath.toRadians(28.), FastMath.toRadians(-30.), 800000.);
346         tc = topoFrame.getTrackingCoordinates(earthSpheric.transform(satPoint), earthSpheric.getBodyFrame(), date);
347 
348         Assertions.assertEquals(FastMath.PI - satPoint.getLongitude(), tc.getAzimuth(), 1.e-5);
349 
350     }
351 
352     @Test
353     void testFieldAzimuthPole() {
354         doTestFieldAzimuthPole(Binary64Field.getInstance());
355     }
356 
357     private <T extends CalculusFieldElement<T>> void doTestFieldAzimuthPole(final Field<T> field) {
358 
359         // zero
360         final T zero = field.getZero();
361 
362         // Surface point at latitude 0
363         final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(89.999), FastMath.toRadians(0.), 0.);
364         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, point, "lon 0 lat 90");
365 
366         // Point at 30 deg longitude
367         // **************************
368         FieldGeodeticPoint<T> satPoint = new FieldGeodeticPoint<>(zero.add(FastMath.toRadians(28.)),
369                                                                   zero.add(FastMath.toRadians(30.)), zero.add(800000.));
370         final FieldAbsoluteDate<T> fieldDate = new FieldAbsoluteDate<>(field, date);
371         FieldTrackingCoordinates<T> tc = topoFrame.getTrackingCoordinates(earthSpheric.transform(satPoint), earthSpheric.getBodyFrame(), fieldDate);
372         Assertions.assertEquals(satPoint.getLongitude().negate().add(FastMath.PI).getReal(),
373                                 tc.getAzimuth().getReal(), 1.e-5);
374 
375         // Point at -30 deg longitude
376         // ***************************
377         satPoint = new FieldGeodeticPoint<>(zero.add(FastMath.toRadians(28.)),
378                                             zero.add(FastMath.toRadians(-30.)), zero.add(800000.));
379         tc = topoFrame.getTrackingCoordinates(earthSpheric.transform(satPoint), earthSpheric.getBodyFrame(), fieldDate);
380 
381         Assertions.assertEquals(satPoint.getLongitude().negate().add(FastMath.PI).getReal(),
382                                 tc.getAzimuth().getReal(), 1.e-5);
383 
384     }
385 
386     @Test
387     void testGetVelocity() {
388         // GIVEN
389         final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(45.), FastMath.toRadians(5.), 0.);
390         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, point, "lon 5 lat 45");
391         final Frame frame = FramesFactory.getGCRF();
392         // WHEN
393         final Vector3D velocity = topoFrame.getVelocity(date, frame);
394         // THEN
395         Assertions.assertEquals(topoFrame.getPVCoordinates(date, frame).getVelocity(), velocity);
396     }
397 
398     @Test
399     void testGetPVCoordinates() {
400         // GIVEN
401         final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(45.), FastMath.toRadians(5.), 0.);
402         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, point, "lon 5 lat 45");
403         final Frame frame = FramesFactory.getGCRF();
404         // WHEN
405         final PVCoordinates pv = topoFrame.getPVCoordinates(date, frame);
406         // THEN
407         final Vector3D position = topoFrame.getPosition(date, frame);
408         Assertions.assertEquals(position, pv.getPosition());
409     }
410 
411     @Test
412     void testFieldGetPVCoordinates() {
413         // GIVEN
414         final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(45.), FastMath.toRadians(5.), 0.);
415         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, point, "lon 5 lat 45");
416         final Frame frame = FramesFactory.getGCRF();
417         final FieldAbsoluteDate<Binary64> fieldDate = new FieldAbsoluteDate<>(Binary64Field.getInstance(), date);
418         // WHEN
419         final FieldPVCoordinates<Binary64> pv = topoFrame.getPVCoordinates(fieldDate, frame);
420         // THEN
421         final FieldVector3D<Binary64> position = topoFrame.getPosition(fieldDate, frame);
422         Assertions.assertEquals(position, pv.getPosition());
423     }
424 
425     @Test
426     void testDoppler() {
427 
428         // Surface point at latitude 45, longitude 5
429         final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(45.), FastMath.toRadians(5.), 0.);
430         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, point, "lon 5 lat 45");
431 
432         // Point at 30 deg longitude
433         // ***************************
434         final CircularOrbit orbit =
435             new CircularOrbit(7178000.0, 0.5e-8, -0.5e-8, FastMath.toRadians(50.), FastMath.toRadians(120.),
436                                    FastMath.toRadians(90.), PositionAngleType.MEAN,
437                                    FramesFactory.getEME2000(), date, mu);
438 
439         // Transform satellite position to position/velocity parameters in body frame
440         final Transform eme2000ToItrf = FramesFactory.getEME2000().getTransformTo(earthSpheric.getBodyFrame(), date);
441         final PVCoordinates pvSatItrf = eme2000ToItrf.transformPVCoordinates(orbit.getPVCoordinates());
442 
443         // Compute range rate directly
444         //********************************************
445         final double dop = topoFrame.getRangeRate(pvSatItrf, earthSpheric.getBodyFrame(), date);
446 
447         // Compare to finite difference computation (2 points)
448         //*****************************************************
449         final double dt = 0.1;
450         KeplerianPropagator extrapolator = new KeplerianPropagator(orbit);
451 
452         // Extrapolate satellite position a short while after reference date
453         AbsoluteDate dateP = date.shiftedBy(dt);
454         Transform j2000ToItrfP = FramesFactory.getEME2000().getTransformTo(earthSpheric.getBodyFrame(), dateP);
455         SpacecraftState orbitP = extrapolator.propagate(dateP);
456         Vector3D satPointGeoP = j2000ToItrfP.transformPosition(orbitP.getPosition());
457 
458         // Retropolate satellite position a short while before reference date
459         AbsoluteDate dateM = date.shiftedBy(-dt);
460         Transform j2000ToItrfM = FramesFactory.getEME2000().getTransformTo(earthSpheric.getBodyFrame(), dateM);
461         SpacecraftState orbitM = extrapolator.propagate(dateM);
462         Vector3D satPointGeoM = j2000ToItrfM.transformPosition(orbitM.getPosition());
463 
464         // Compute ranges at both instants
465         double rangeP = topoFrame.getTrackingCoordinates(satPointGeoP, earthSpheric.getBodyFrame(), dateP).getRange();
466         double rangeM = topoFrame.getTrackingCoordinates(satPointGeoM, earthSpheric.getBodyFrame(), dateM).getRange();
467         final double dopRef2 = (rangeP - rangeM) / (2. * dt);
468         Assertions.assertEquals(dopRef2, dop, 1.e-3);
469 
470     }
471 
472     @Test
473     void testFieldDoppler() {
474         doTestFieldDoppler(Binary64Field.getInstance());
475     }
476 
477     private <T extends CalculusFieldElement<T>> void doTestFieldDoppler(final Field<T> field) {
478 
479         // zero
480         final T zero = field.getZero();
481 
482         // Surface point at latitude 45, longitude 5
483         final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(45.), FastMath.toRadians(5.), 0.);
484         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, point, "lon 5 lat 45");
485 
486         // Field date
487         final FieldAbsoluteDate<T> fieldDate = new FieldAbsoluteDate<>(field, date);
488 
489         // Point at 30 deg longitude
490         // ***************************
491         final FieldCircularOrbit<T> orbit =
492             new FieldCircularOrbit<>(zero.add(7178000.0), zero.add(0.5e-8), zero.add(-0.5e-8),
493                                      zero.add(FastMath.toRadians(50.)), zero.add(FastMath.toRadians(120.)),
494                                      zero.add(FastMath.toRadians(90.)), PositionAngleType.MEAN,
495                                      FramesFactory.getEME2000(), fieldDate, zero.add(mu));
496 
497         // Transform satellite position to position/velocity parameters in body frame
498         final FieldTransform<T> eme2000ToItrf = FramesFactory.getEME2000().getTransformTo(earthSpheric.getBodyFrame(), fieldDate);
499         final FieldPVCoordinates<T> pvSatItrf = eme2000ToItrf.transformPVCoordinates(orbit.getPVCoordinates());
500 
501         // Compute range rate directly
502         //********************************************
503         final T dop = topoFrame.getRangeRate(pvSatItrf, earthSpheric.getBodyFrame(), fieldDate);
504 
505         // Compare to finite difference computation (2 points)
506         //*****************************************************
507         final double dt = 0.1;
508         FieldKeplerianPropagator<T> extrapolator = new FieldKeplerianPropagator<>(orbit);
509 
510         // Extrapolate satellite position a short while after reference date
511         FieldAbsoluteDate<T> dateP = fieldDate.shiftedBy(dt);
512         FieldTransform<T> j2000ToItrfP = FramesFactory.getEME2000().getTransformTo(earthSpheric.getBodyFrame(), dateP);
513         FieldSpacecraftState<T> orbitP = extrapolator.propagate(dateP);
514         FieldVector3D<T> satPointGeoP = j2000ToItrfP.transformPosition(orbitP.getPosition());
515 
516         // Retropolate satellite position a short while before reference date
517         FieldAbsoluteDate<T> dateM = fieldDate.shiftedBy(-dt);
518         FieldTransform<T> j2000ToItrfM = FramesFactory.getEME2000().getTransformTo(earthSpheric.getBodyFrame(), dateM);
519         FieldSpacecraftState<T> orbitM = extrapolator.propagate(dateM);
520         FieldVector3D<T> satPointGeoM = j2000ToItrfM.transformPosition(orbitM.getPosition());
521 
522         // Compute ranges at both instants
523         T rangeP = topoFrame.getTrackingCoordinates(satPointGeoP, earthSpheric.getBodyFrame(), dateP).getRange();
524         T rangeM = topoFrame.getTrackingCoordinates(satPointGeoM, earthSpheric.getBodyFrame(), dateM).getRange();
525         final T dopRef2 = (rangeP.subtract(rangeM)).divide(2. * dt);
526         Assertions.assertEquals(dopRef2.getReal(), dop.getReal(), 1.e-3);
527 
528     }
529 
530     @Test
531     void testEllipticEarth()  {
532 
533         // Elliptic earth shape
534         final OneAxisEllipsoid earthElliptic =
535             new OneAxisEllipsoid(6378136.460, 1 / 298.257222101, itrf);
536 
537         // Satellite point
538         // Caution !!! Sat point target shall be the same whatever earth shape chosen !!
539         final GeodeticPoint satPointGeo = new GeodeticPoint(FastMath.toRadians(30.), FastMath.toRadians(15.), 800000.);
540         final Vector3D satPoint = earthElliptic.transform(satPointGeo);
541 
542         // ****************************
543         // Test at equatorial position
544         // ****************************
545         GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(0.), FastMath.toRadians(5.), 0.);
546         TopocentricFrame topoElliptic  = new TopocentricFrame(earthElliptic, point, "elliptic, equatorial lon 5");
547         TopocentricFrame topoSpheric   = new TopocentricFrame(earthSpheric, point, "spheric, equatorial lon 5");
548         TrackingCoordinates tcElliptic = topoElliptic.getTrackingCoordinates(satPoint, earthElliptic.getBodyFrame(), date);
549         TrackingCoordinates tcSpheric  = topoSpheric.getTrackingCoordinates(satPoint, earthElliptic.getBodyFrame(), date);
550 
551         // Compare azimuth/elevation/range of satellite point : shall be strictly identical
552         // ***************************************************
553         Assertions.assertEquals(tcElliptic.getAzimuth(),   tcSpheric.getAzimuth(),   Utils.epsilonAngle);
554         Assertions.assertEquals(tcElliptic.getElevation(), tcSpheric.getElevation(), Utils.epsilonAngle);
555         Assertions.assertEquals(tcElliptic.getRange(),     tcSpheric.getRange(),     Utils.epsilonTest);
556 
557         // Infinite point separated by -20 deg in longitude
558         // *************************************************
559         GeodeticPoint infPointGeo = new GeodeticPoint(FastMath.toRadians(0.), FastMath.toRadians(-15.), 1000000000.);
560         Vector3D infPoint = earthElliptic.transform(infPointGeo);
561         tcElliptic = topoElliptic.getTrackingCoordinates(infPoint, earthElliptic.getBodyFrame(), date);
562 
563         // Azimuth = pi/2
564         Assertions.assertEquals(3*FastMath.PI/2., tcElliptic.getAzimuth(), Utils.epsilonAngle);
565 
566         // Elevation = pi/2 - longitude difference
567         Assertions.assertEquals(FastMath.PI/2. - FastMath.abs(point.getLongitude() - infPointGeo.getLongitude()),
568                                 tcElliptic.getElevation(), 1.e-2);
569 
570         // Infinite point separated by +20 deg in longitude
571         // *************************************************
572         infPointGeo = new GeodeticPoint(FastMath.toRadians(0.), FastMath.toRadians(25.), 1000000000.);
573         infPoint = earthElliptic.transform(infPointGeo);
574         tcElliptic = topoElliptic.getTrackingCoordinates(infPoint, earthElliptic.getBodyFrame(), date);
575 
576         // Azimuth = pi/2
577         Assertions.assertEquals(FastMath.PI/2., tcElliptic.getAzimuth(), Utils.epsilonAngle);
578 
579         // Site = pi/2 - longitude difference
580         Assertions.assertEquals(FastMath.PI/2. - FastMath.abs(point.getLongitude() - infPointGeo.getLongitude()),
581                                 tcElliptic.getElevation(), 1.e-2);
582 
583         // ************************
584         // Test at polar position
585         // ************************
586         point = new GeodeticPoint(FastMath.toRadians(89.999), FastMath.toRadians(0.), 0.);
587         topoSpheric  = new TopocentricFrame(earthSpheric, point, "lon 0 lat 90");
588         topoElliptic = new TopocentricFrame(earthElliptic, point, "lon 0 lat 90");
589 
590         // Compare azimuth/elevation/range of satellite point : slight difference due to earth flatness
591         // ***************************************************
592         tcElliptic = topoElliptic.getTrackingCoordinates(satPoint, earthElliptic.getBodyFrame(), date);
593         tcSpheric  = topoSpheric.getTrackingCoordinates(satPoint, earthSpheric.getBodyFrame(), date);
594         Assertions.assertEquals(tcElliptic.getAzimuth(),   tcSpheric.getAzimuth(),   1.e-7);
595         Assertions.assertEquals(tcElliptic.getElevation(), tcSpheric.getElevation(), 1.e-2);
596         Assertions.assertEquals(tcElliptic.getRange(),     tcSpheric.getRange(),     20.e+3);
597 
598 
599         // *********************
600         // Test at any position
601         // *********************
602         point = new GeodeticPoint(FastMath.toRadians(60), FastMath.toRadians(30.), 0.);
603         topoSpheric  = new TopocentricFrame(earthSpheric, point, "lon 10 lat 45");
604         topoElliptic = new TopocentricFrame(earthElliptic, point, "lon 10 lat 45");
605 
606         // Compare azimuth/elevation/range of satellite point : slight difference
607         // ***************************************************
608         tcElliptic = topoElliptic.getTrackingCoordinates(satPoint, earthElliptic.getBodyFrame(), date);
609         tcSpheric  = topoSpheric.getTrackingCoordinates(satPoint, earthSpheric.getBodyFrame(), date);
610         Assertions.assertEquals(tcElliptic.getAzimuth(),   tcSpheric.getAzimuth(),   1.e-2);
611         Assertions.assertEquals(tcElliptic.getElevation(), tcSpheric.getElevation(), 1.e-2);
612         Assertions.assertEquals(tcElliptic.getRange(),     tcSpheric.getRange(),     20.e+3);
613 
614     }
615 
616     @Test
617     void testFieldEllipticEarth() {
618         doTestFieldEllipticEarth(Binary64Field.getInstance());
619     }
620 
621     private <T extends CalculusFieldElement<T>> void doTestFieldEllipticEarth(final Field<T> field)  {
622 
623         // zero
624         final T zero = field.getZero();
625 
626         // Elliptic earth shape
627         final OneAxisEllipsoid earthElliptic =
628             new OneAxisEllipsoid(6378136.460, 1 / 298.257222101, itrf);
629 
630         // Field date
631         final FieldAbsoluteDate<T> fieldDate = new FieldAbsoluteDate<>(field, date);
632 
633         // Satellite point
634         // Caution !!! Sat point target shall be the same whatever earth shape chosen !!
635         final FieldGeodeticPoint<T> satPointGeo = new FieldGeodeticPoint<>(zero.add(FastMath.toRadians(30.)),
636                                                                            zero.add(FastMath.toRadians(15.)),
637                                                                            zero.add(800000.));
638         final FieldVector3D<T> satPoint = earthElliptic.transform(satPointGeo);
639 
640         // ****************************
641         // Test at equatorial position
642         // ****************************
643         GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(0.), FastMath.toRadians(5.), 0.);
644         TopocentricFrame topoElliptic = new TopocentricFrame(earthElliptic, point, "elliptic, equatorial lon 5");
645         TopocentricFrame topoSpheric = new TopocentricFrame(earthSpheric, point, "spheric, equatorial lon 5");
646 
647         // Compare azimuth/elevation/range of satellite point : shall be strictly identical
648         // ***************************************************
649         FieldTrackingCoordinates<T> tcElliptic =
650                         topoElliptic.getTrackingCoordinates(satPoint, earthElliptic.getBodyFrame(), fieldDate);
651         FieldTrackingCoordinates<T> tcSpheric =
652                         topoSpheric.getTrackingCoordinates(satPoint, earthSpheric.getBodyFrame(), fieldDate);
653         Assertions.assertEquals(tcElliptic.getAzimuth().getReal(),   tcSpheric.getAzimuth().getReal(),   Utils.epsilonAngle);
654         Assertions.assertEquals(tcElliptic.getElevation().getReal(), tcSpheric.getElevation().getReal(), Utils.epsilonAngle);
655         Assertions.assertEquals(tcElliptic.getRange().getReal(),     tcSpheric.getRange().getReal(),     Utils.epsilonTest);
656 
657         // Infinite point separated by -20 deg in longitude
658         // *************************************************
659         FieldGeodeticPoint<T> infPointGeo = new FieldGeodeticPoint<>(zero.add(FastMath.toRadians(0.)),
660                                                                      zero.add(FastMath.toRadians(-15.)),
661                                                                      zero.add(1000000000.));
662         FieldVector3D<T> infPoint = earthElliptic.transform(infPointGeo);
663         tcElliptic = topoElliptic.getTrackingCoordinates(infPoint, earthElliptic.getBodyFrame(), fieldDate);
664 
665         // Azimuth = pi/2
666         Assertions.assertEquals(3*FastMath.PI/2., tcElliptic.getAzimuth().getReal(), Utils.epsilonAngle);
667 
668         // Elevation = pi/2 - longitude difference
669         Assertions.assertEquals(FastMath.abs(infPointGeo.getLongitude().negate().add(point.getLongitude())).negate().add(FastMath.PI/2.).getReal(),
670                                 tcElliptic.getElevation().getReal(),
671                                 1.e-2);
672 
673         // Infinite point separated by +20 deg in longitude
674         // *************************************************
675         infPointGeo = new FieldGeodeticPoint<>(zero.add(FastMath.toRadians(0.)),
676                                                zero.add(FastMath.toRadians(25.)),
677                                                zero.add(1000000000.));
678         infPoint = earthElliptic.transform(infPointGeo);
679         tcElliptic = topoElliptic.getTrackingCoordinates(infPoint, earthElliptic.getBodyFrame(), fieldDate);
680 
681         // Azimuth = pi/2
682         Assertions.assertEquals(FastMath.PI/2., tcElliptic.getAzimuth().getReal(), Utils.epsilonAngle);
683 
684         // Elevation = pi/2 - longitude difference
685         Assertions.assertEquals(FastMath.abs(infPointGeo.getLongitude().negate().add(point.getLongitude())).negate().add(FastMath.PI/2.).getReal(),
686                                 tcElliptic.getElevation().getReal(),
687                                 1.e-2);
688 
689         // ************************
690         // Test at polar position
691         // ************************
692         point = new GeodeticPoint(FastMath.toRadians(89.999), FastMath.toRadians(0.), 0.);
693         topoSpheric  = new TopocentricFrame(earthSpheric, point, "lon 0 lat 90");
694         topoElliptic = new TopocentricFrame(earthElliptic, point, "lon 0 lat 90");
695         tcElliptic = topoElliptic.getTrackingCoordinates(satPoint, earthElliptic.getBodyFrame(), fieldDate);
696         tcSpheric = topoSpheric.getTrackingCoordinates(satPoint, earthSpheric.getBodyFrame(), fieldDate);
697 
698         // Compare azimuth/elevation/range of satellite point : slight difference due to earth flatness
699         // ***************************************************
700         Assertions.assertEquals(tcElliptic.getAzimuth().getReal(),   tcSpheric.getAzimuth().getReal(),   1.e-7);
701         Assertions.assertEquals(tcElliptic.getElevation().getReal(), tcSpheric.getElevation().getReal(), 1.e-2);
702         Assertions.assertEquals(tcElliptic.getRange().getReal(),     tcSpheric.getRange().getReal(),     20.e+3);
703 
704 
705         // *********************
706         // Test at any position
707         // *********************
708         point = new GeodeticPoint(FastMath.toRadians(60), FastMath.toRadians(30.), 0.);
709         topoSpheric  = new TopocentricFrame(earthSpheric, point, "lon 10 lat 45");
710         topoElliptic = new TopocentricFrame(earthElliptic, point, "lon 10 lat 45");
711         tcElliptic = topoElliptic.getTrackingCoordinates(satPoint, earthElliptic.getBodyFrame(), fieldDate);
712         tcSpheric = topoSpheric.getTrackingCoordinates(satPoint, earthSpheric.getBodyFrame(), fieldDate);
713 
714         // Compare azimuth/elevation/range of satellite point : slight difference
715         // ***************************************************
716         Assertions.assertEquals(tcElliptic.getAzimuth().getReal(),   tcSpheric.getAzimuth().getReal(),   1.e-2);
717         Assertions.assertEquals(tcElliptic.getElevation().getReal(), tcSpheric.getElevation().getReal(), 1.e-2);
718         Assertions.assertEquals(tcElliptic.getRange().getReal(),     tcSpheric.getRange().getReal(),     20.e+3);
719 
720     }
721 
722     @Test
723     void testPointAtDistance() {
724 
725         RandomGenerator random = new Well1024a(0xa1e6bd5cd0578779l);
726         final OneAxisEllipsoid earth =
727             new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING,
728                                  itrf);
729         final AbsoluteDate date = AbsoluteDate.J2000_EPOCH;
730 
731         for (int i = 0; i < 20; ++i) {
732             // we don't need uniform point on the sphere, just a few different test configurations
733             double latitude  = FastMath.PI * (0.5 - random.nextDouble());
734             double longitude = 2 * FastMath.PI * random.nextDouble();
735             TopocentricFrame topo = new TopocentricFrame(earth,
736                                                          new GeodeticPoint(latitude, longitude, 0.0),
737                                                          "topo");
738             Transform transform = earth.getBodyFrame().getTransformTo(topo, date);
739             for (int j = 0; j < 20; ++j) {
740                 double elevation      = FastMath.PI * (0.5 - random.nextDouble());
741                 double azimuth        = 2 * FastMath.PI * random.nextDouble();
742                 double range          = 500000.0 * (1.0 + random.nextDouble());
743                 Vector3D absolutePoint = earth.transform(topo.pointAtDistance(azimuth, elevation, range));
744                 Vector3D relativePoint = transform.transformPosition(absolutePoint);
745                 TrackingCoordinates rebuiltTC = topo.getTrackingCoordinates(relativePoint, topo, AbsoluteDate.J2000_EPOCH);
746                 Assertions.assertEquals(elevation, rebuiltTC.getElevation(), 1.0e-12);
747                 Assertions.assertEquals(azimuth, MathUtils.normalizeAngle(rebuiltTC.getAzimuth(), azimuth), 1.0e-12);
748                 Assertions.assertEquals(range, rebuiltTC.getRange(), 1.0e-12 * range);
749             }
750         }
751     }
752 
753     @Test
754     void testIssue145() {
755         Frame itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
756         BodyShape earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
757                                                Constants.WGS84_EARTH_FLATTENING,
758                                                itrf);
759         TopocentricFrame staFrame = new TopocentricFrame(earth, new GeodeticPoint(0.0, 0.0, 0.0), "test");
760         GeodeticPoint gp = staFrame.computeLimitVisibilityPoint(Constants.WGS84_EARTH_EQUATORIAL_RADIUS+600000,
761                                                                 0.0, FastMath.toRadians(5.0));
762         Assertions.assertEquals(0.0, gp.getLongitude(), 1.0e-15);
763         Assertions.assertTrue(gp.getLatitude() > 0);
764         Assertions.assertEquals(0.0, staFrame.getNorth().distance(Vector3D.PLUS_K), 1.0e-15);
765 
766     }
767 
768     @Test
769     void testVisibilityCircle() throws IOException {
770 
771         // a few random from International Laser Ranging Service
772         final BodyShape earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
773                                                      Constants.WGS84_EARTH_FLATTENING,
774                                                      FramesFactory.getITRF(IERSConventions.IERS_2010, true));
775         final TopocentricFrame[] ilrs = {
776             new TopocentricFrame(earth,
777                                  new GeodeticPoint(FastMath.toRadians(52.3800), FastMath.toRadians(3.0649), 133.745),
778                                  "Potsdam"),
779             new TopocentricFrame(earth,
780                                  new GeodeticPoint(FastMath.toRadians(36.46273), FastMath.toRadians(-6.20619), 64.0),
781                                  "San Fernando"),
782             new TopocentricFrame(earth,
783                                  new GeodeticPoint(FastMath.toRadians(35.5331), FastMath.toRadians(24.0705), 157.0),
784                                  "Chania")
785         };
786 
787         PolynomialFunction distanceModel =
788                 new PolynomialFunction(new double[] { 7.0892e+05, 3.1913, -8.2181e-07, 1.4033e-13 });
789         for (TopocentricFrame station : ilrs) {
790             for (double altitude = 500000; altitude < 2000000; altitude += 100000) {
791                 for (double azimuth = 0; azimuth < 2 * FastMath.PI; azimuth += 0.05) {
792                     GeodeticPoint p = station.computeLimitVisibilityPoint(Constants.WGS84_EARTH_EQUATORIAL_RADIUS + altitude,
793                                                                           azimuth, FastMath.toRadians(5.0));
794                     double d = station.getTrackingCoordinates(earth.transform(p), earth.getBodyFrame(), AbsoluteDate.J2000_EPOCH).getRange();
795                     Assertions.assertEquals(distanceModel.value(altitude), d, 40000.0);
796                 }
797             }
798         }
799 
800     }
801 
802     private static Stream<Arguments> testGetTopocentricCoordinatesValues() {
803         return Stream.of(
804                 Arguments.of(0,0,1),
805                 Arguments.of(0, FastMath.PI / 2, 1),
806                 Arguments.of(FastMath.PI, FastMath.PI / 3, 10),
807                 Arguments.of(3 * FastMath.PI / 2, FastMath.PI / 4, 1000),
808                 Arguments.of(FastMath.PI / 2, -FastMath.PI / 6, 500),
809                 Arguments.of(FastMath.PI / 7, -FastMath.PI / 5, 100)
810         );
811     }
812 
813     @ParameterizedTest
814     @MethodSource("testGetTopocentricCoordinatesValues")
815     void testGetTopocentricCoordinates(final double az, final double el, final double r) {
816         final TrackingCoordinates coords = new TrackingCoordinates(az,el,r);
817         final Vector3D topoPos = TopocentricFrame.getTopocentricPosition(coords);
818 
819         final SinCos sinCosAz = FastMath.sinCos(az);
820         final SinCos sinCosEl = FastMath.sinCos(el);
821         final double expectedX, expectedY, expectedZ;
822         expectedX = sinCosAz.sin() * sinCosEl.cos() * r;
823         expectedY = sinCosAz.cos() * sinCosEl.cos() * r;
824         expectedZ = sinCosEl.sin() * r;
825 
826         final double distanceTolerance = 1e-7;
827         Assertions.assertEquals(expectedX,topoPos.getX(), distanceTolerance);
828         Assertions.assertEquals(expectedY,topoPos.getY(), distanceTolerance);
829         Assertions.assertEquals(expectedZ,topoPos.getZ(), distanceTolerance);
830     }
831 
832     @ParameterizedTest
833     @MethodSource("testGetTopocentricCoordinatesValues")
834     void testGetFieldTopocentricCoordinates(final double azimuth, final double elevation, final double range) {
835         final Binary64 az, el, r;
836         az = new Binary64(azimuth);
837         el = new Binary64(elevation);
838         r = new Binary64(range);
839 
840         final FieldTrackingCoordinates<Binary64> coords = new FieldTrackingCoordinates<>(az, el, r);
841         final FieldVector3D<Binary64> topoPos = TopocentricFrame.getTopocentricPosition(coords);
842 
843         final FieldSinCos<Binary64> sinCosAz = FastMath.sinCos(az);
844         final FieldSinCos<Binary64> sinCosEl = FastMath.sinCos(el);
845         final Binary64 expectedX, expectedY, expectedZ;
846         expectedX = r.multiply(sinCosAz.sin().multiply(sinCosEl.cos()));
847         expectedY = r.multiply(sinCosEl.cos().multiply(sinCosAz.cos()));
848         expectedZ = r.multiply(sinCosEl.sin());
849 
850         final double distanceTolerance = 1e-7;
851         Assertions.assertEquals(expectedX.getReal(), topoPos.getX().getReal(), distanceTolerance);
852         Assertions.assertEquals(expectedY.getReal(), topoPos.getY().getReal(), distanceTolerance);
853         Assertions.assertEquals(expectedZ.getReal(), topoPos.getZ().getReal(), distanceTolerance);
854     }
855 
856     @ParameterizedTest
857     @MethodSource("testGetTopocentricCoordinatesValues")
858     void testInverseGetTopocentricCoordinates(double az, double el, double r) {
859         final TrackingCoordinates expectedCoords = new TrackingCoordinates(az, el, r);
860 
861         final Vector3D point = TopocentricFrame.getTopocentricPosition(expectedCoords);
862         final GeodeticPoint geodeticPoint = new GeodeticPoint(0., 0., 0.);
863         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, geodeticPoint, "geodeticPoint");
864         final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
865         final TrackingCoordinates coords = topoFrame.getTrackingCoordinates(point, topoFrame, date);
866 
867         final double angularTolerance = 1e-12;
868         Assertions.assertEquals(expectedCoords.getAzimuth(), coords.getAzimuth(), angularTolerance);
869         Assertions.assertEquals(expectedCoords.getElevation(), coords.getElevation(), angularTolerance);
870         Assertions.assertEquals(expectedCoords.getRange(), coords.getRange());
871     }
872 
873     @Test
874     void testGetTrackingCoordinates() {
875         // GIVEN
876         final GeodeticPoint geodeticPoint = new GeodeticPoint(0., 0., 0.);
877         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, geodeticPoint, "geodeticPoint");
878         final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
879         final Frame frame = FramesFactory.getGCRF();
880         final Vector3D point = new Vector3D(1., 2., -3.);
881         // WHEN
882         final TrackingCoordinates trackingCoordinates = topoFrame.getTrackingCoordinates(point, frame, date);
883         // THEN
884         final double expectedElevation = topoFrame.getElevation(point, frame, date);
885         final double expectedAzimuth = topoFrame.getAzimuth(point, frame, date);
886         final double expectedRange = topoFrame.getRange(point, frame, date);
887         Assertions.assertEquals(expectedElevation, trackingCoordinates.getElevation());
888         Assertions.assertEquals(expectedAzimuth, trackingCoordinates.getAzimuth());
889         Assertions.assertEquals(expectedRange, trackingCoordinates.getRange());
890     }
891 
892     @Test
893     void testGetFieldTrackingCoordinates() {
894         // GIVEN
895         final GeodeticPoint geodeticPoint = new GeodeticPoint(0., 0., 0.);
896         final TopocentricFrame topoFrame = new TopocentricFrame(earthSpheric, geodeticPoint, "geodeticPoint");
897         final ComplexField field = ComplexField.getInstance();
898         final FieldAbsoluteDate<Complex> date = FieldAbsoluteDate.getArbitraryEpoch(field);
899         final Frame frame = FramesFactory.getGCRF();
900         final FieldVector3D<Complex> point = new FieldVector3D<>(field, new Vector3D(1., 2., -3.));
901         // WHEN
902         final FieldTrackingCoordinates<Complex> trackingCoordinates = topoFrame.getTrackingCoordinates(point, frame, date);
903         // THEN
904         final Complex expectedElevation = topoFrame.getElevation(point, frame, date);
905         final Complex expectedAzimuth = topoFrame.getAzimuth(point, frame, date);
906         final Complex expectedRange = topoFrame.getRange(point, frame, date);
907         Assertions.assertEquals(expectedElevation, trackingCoordinates.getElevation());
908         Assertions.assertEquals(expectedAzimuth, trackingCoordinates.getAzimuth());
909         Assertions.assertEquals(expectedRange, trackingCoordinates.getRange());
910     }
911 
912     @BeforeEach
913     public void setUp() {
914 
915         Utils.setDataRoot("regular-data");
916 
917         // Reference frame = ITRF
918         itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
919 
920         // Elliptic earth shape
921         earthSpheric = new OneAxisEllipsoid(6378136.460, 0., itrf);
922 
923         // Reference date
924         date = new AbsoluteDate(new DateComponents(2008, 04, 07),
925                                 TimeComponents.H00,
926                                 TimeScalesFactory.getUTC());
927 
928         // Body mu
929         mu = 3.9860047e14;
930 
931     }
932 
933     @AfterEach
934     public void tearDown() {
935         date = null;
936         itrf = null;
937         earthSpheric = null;
938     }
939 
940 
941 }