1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17 package org.orekit.rugged.intersection;
18
19
20 import java.io.File;
21 import java.net.URISyntaxException;
22
23 import org.hipparchus.geometry.euclidean.threed.Line;
24 import org.hipparchus.geometry.euclidean.threed.Rotation;
25 import org.hipparchus.geometry.euclidean.threed.Vector3D;
26 import org.hipparchus.util.FastMath;
27 import org.junit.jupiter.api.AfterEach;
28 import org.junit.jupiter.api.Assertions;
29 import org.junit.jupiter.api.BeforeEach;
30 import org.junit.jupiter.api.Test;
31 import org.orekit.attitudes.Attitude;
32 import org.orekit.bodies.GeodeticPoint;
33 import org.orekit.data.DataContext;
34 import org.orekit.data.DirectoryCrawler;
35 import org.orekit.frames.FramesFactory;
36 import org.orekit.frames.Transform;
37 import org.orekit.orbits.CartesianOrbit;
38 import org.orekit.propagation.SpacecraftState;
39 import org.orekit.rugged.intersection.duvenhage.MinMaxTreeTile;
40 import org.orekit.rugged.intersection.duvenhage.MinMaxTreeTileFactory;
41 import org.orekit.rugged.raster.CliffsElevationUpdater;
42 import org.orekit.rugged.raster.TileUpdater;
43 import org.orekit.rugged.raster.VolcanicConeElevationUpdater;
44 import org.orekit.rugged.utils.ExtendedEllipsoid;
45 import org.orekit.time.AbsoluteDate;
46 import org.orekit.time.TimeScalesFactory;
47 import org.orekit.utils.Constants;
48 import org.orekit.utils.IERSConventions;
49 import org.orekit.utils.PVCoordinates;
50
51 public abstract class AbstractAlgorithmTest {
52
53 protected abstract IntersectionAlgorithm createAlgorithm(TileUpdater updater, int maxCachedTiles, boolean isOverlappingTiles);
54
55 @Test
56 public void testMayonVolcanoOnSubTileCorner() {
57
58 setUpMayonVolcanoContext();
59
60
61
62
63 final double latitude = FastMath.toRadians(13.27);
64 final double longitude = FastMath.toRadians(123.68);
65 MinMaxTreeTile tile = new MinMaxTreeTileFactory().createTile();
66 updater.updateTile(latitude, longitude, tile);
67 double altitude = tile.interpolateElevation(latitude, longitude);
68 final GeodeticPoint groundGP = new GeodeticPoint(latitude, longitude, altitude);
69 Vector3D groundP = earth.transform(groundGP);
70
71 final IntersectionAlgorithm algorithm = createAlgorithm(updater, 8, true);
72
73
74 Transform earthToSpacecraft = new Transform(state.getDate(),
75 earth.getBodyFrame().getTransformTo(state.getFrame(), state.getDate()),
76 state.toTransform());
77 Vector3D pointInSpacecraftFrame = earthToSpacecraft.transformPosition(groundP);
78 Assertions.assertEquals( 0.000, pointInSpacecraftFrame.getX(), 1.0e-3);
79 Assertions.assertEquals(-87754.914, pointInSpacecraftFrame.getY(), 1.0e-3);
80 Assertions.assertEquals(790330.254, pointInSpacecraftFrame.getZ(), 1.0e-3);
81
82
83 Vector3D position = state.getPVCoordinates(earth.getBodyFrame()).getPosition();
84 Vector3D los = groundP.subtract(position);
85 GeodeticPoint result = algorithm.refineIntersection(earth, position, los,
86 algorithm.intersection(earth, position, los));
87 checkIntersection(position, los, result);
88 Assertions.assertEquals(0.0, groundP.distance(earth.transform(result)), 2.0e-9);
89
90 }
91
92 @Test
93 public void testMayonVolcanoWithinPixel() {
94
95 setUpMayonVolcanoContext();
96
97 final double latitude = FastMath.toRadians(13.2696);
98 final double longitude = FastMath.toRadians(123.6803);
99 MinMaxTreeTile tile = new MinMaxTreeTileFactory().createTile();
100 updater.updateTile(latitude, longitude, tile);
101 double altitude = tile.interpolateElevation(latitude, longitude);
102 final GeodeticPoint groundGP = new GeodeticPoint(latitude, longitude, altitude);
103 Vector3D groundP = earth.transform(groundGP);
104
105 final IntersectionAlgorithm algorithm = createAlgorithm(updater, 8, true);
106
107
108 Vector3D position = state.getPVCoordinates(earth.getBodyFrame()).getPosition();
109 Vector3D los = groundP.subtract(position);
110 GeodeticPoint result = algorithm.refineIntersection(earth, position, los,
111 algorithm.intersection(earth, position, los));
112 checkIntersection(position, los, result);
113 Assertions.assertEquals(0.0, groundP.distance(earth.transform(result)), 2.0e-9);
114
115 }
116
117 @Test
118 public void testCliffsOfMoher()
119 {
120
121 setUpCliffsOfMoherContext();
122
123
124 final double latitude = FastMath.toRadians(52.98045);
125 final double longitude = FastMath.toRadians(-9.421826199814143);
126 MinMaxTreeTile tile = new MinMaxTreeTileFactory().createTile();
127 updater.updateTile(latitude, longitude, tile);
128 double altitude = tile.interpolateElevation(latitude, longitude);
129 final GeodeticPoint groundGP = new GeodeticPoint(latitude, longitude, altitude);
130 Vector3D groundP = earth.transform(groundGP);
131
132 final IntersectionAlgorithm algorithm = createAlgorithm(updater, 8, true);
133 Assertions.assertEquals( 0.0, algorithm.getElevation(latitude, longitude - 2.0e-5), 1.0e-6);
134 Assertions.assertEquals(120.0, algorithm.getElevation(latitude, longitude + 2.0e-5), 1.0e-6);
135
136
137 Transform earthToSpacecraft = new Transform(state.getDate(),
138 earth.getBodyFrame().getTransformTo(state.getFrame(), state.getDate()),
139 state.toTransform());
140 Vector3D pointInSpacecraftFrame = earthToSpacecraft.transformPosition(groundP);
141 Assertions.assertEquals( 0.000, pointInSpacecraftFrame.getX(), 1.0e-3);
142 Assertions.assertEquals( 66702.419, pointInSpacecraftFrame.getY(), 1.0e-3);
143 Assertions.assertEquals(796873.178, pointInSpacecraftFrame.getZ(), 1.0e-3);
144
145 Vector3D position = state.getPVCoordinates(earth.getBodyFrame()).getPosition();
146 Vector3D los = groundP.subtract(position);
147 GeodeticPoint result = algorithm.refineIntersection(earth, position, los,
148 algorithm.intersection(earth, position, los));
149 checkIntersection(position, los, result);
150 Assertions.assertEquals(0.0, groundP.distance(earth.transform(result)), 2.0e-9);
151
152 }
153
154 protected void checkIntersection(Vector3D position, Vector3D los, GeodeticPoint intersection) {
155
156
157 Line line = new Line(position, new Vector3D(1, position, 1e6, los), 1.0e-12);
158 Assertions.assertEquals(0.0, line.distance(earth.transform(intersection)), 3.0e-9);
159
160
161 MinMaxTreeTile tile = new MinMaxTreeTileFactory().createTile();
162 updater.updateTile(intersection.getLatitude(), intersection.getLongitude(), tile);
163 Assertions.assertEquals(tile.interpolateElevation(intersection.getLatitude(), intersection.getLongitude()),
164 intersection.getAltitude(), 2.0e-9);
165
166 }
167
168 protected void setUpMayonVolcanoContext()
169 {
170
171
172 GeodeticPoint summit =
173 new GeodeticPoint(FastMath.toRadians(13.25667), FastMath.toRadians(123.685), 2463.0);
174 updater = new VolcanicConeElevationUpdater(summit,
175 FastMath.toRadians(30.0), 16.0,
176 FastMath.toRadians(1.0), 1201);
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192 AbsoluteDate crossing = new AbsoluteDate("2012-01-06T02:27:15.942757185", TimeScalesFactory.getUTC());
193 state = new SpacecraftState(new CartesianOrbit(new PVCoordinates(new Vector3D( -649500.423763743,
194 -6943715.537565755,
195 1657929.137063380),
196 new Vector3D(-1305.453711368668,
197 -1600.627551928136,
198 -7167.286855869801)),
199 FramesFactory.getEME2000(),
200 crossing,
201 Constants.EIGEN5C_EARTH_MU),
202 new Attitude(crossing,
203 FramesFactory.getEME2000(),
204 new Rotation(-0.40904880353552850,
205 0.46125295378582530,
206 -0.63525007056319790,
207 -0.46516893361386025,
208 true),
209 new Vector3D(-7.048568391860185e-05,
210 -1.043582650222194e-03,
211 1.700400341147713e-05),
212 Vector3D.ZERO));
213
214 }
215
216 protected void setUpCliffsOfMoherContext()
217 {
218
219
220 GeodeticPoint north = new GeodeticPoint(FastMath.toRadians(52.9984),
221 FastMath.toRadians(-9.4072),
222 120.0);
223 GeodeticPoint south = new GeodeticPoint(FastMath.toRadians(52.9625),
224 FastMath.toRadians(-9.4369),
225 120.0);
226
227
228 updater = new CliffsElevationUpdater(north, south,
229 120.0, 0.0,
230 FastMath.toRadians(0.015), 101);
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246 AbsoluteDate crossing = new AbsoluteDate("2012-01-07T11:50:04.935272115", TimeScalesFactory.getUTC());
247 state = new SpacecraftState(new CartesianOrbit(new PVCoordinates(new Vector3D( 412324.544397459,
248 -4325872.329311633,
249 5692124.593989491),
250 new Vector3D(-1293.174701214779,
251 -5900.764863603793,
252 -4378.671036383179)),
253 FramesFactory.getEME2000(),
254 crossing,
255 Constants.EIGEN5C_EARTH_MU),
256 new Attitude(crossing,
257 FramesFactory.getEME2000(),
258 new Rotation(-0.17806699079182878,
259 0.60143347387211290,
260 -0.73251248177468900,
261 -0.26456641385623986,
262 true),
263 new Vector3D(-4.289600857433520e-05,
264 -1.039151496480297e-03,
265 5.811423736843181e-05),
266 Vector3D.ZERO));
267
268 }
269
270 @BeforeEach
271 public void setUp() throws URISyntaxException {
272
273 String path = getClass().getClassLoader().getResource("orekit-data").toURI().getPath();
274 DataContext.getDefault().getDataProvidersManager().addProvider(new DirectoryCrawler(new File(path)));
275 earth = new ExtendedEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
276 Constants.WGS84_EARTH_FLATTENING,
277 FramesFactory.getITRF(IERSConventions.IERS_2010, true));
278 }
279
280 @AfterEach
281 public void tearDown() {
282 earth = null;
283 updater = null;
284 state = null;
285 }
286
287 public ExtendedEllipsoid earth;
288 protected TileUpdater updater;
289 protected SpacecraftState state;
290
291 }