1   /* Copyright 2013-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.rugged.intersection.duvenhage;
18  
19  
20  import org.hipparchus.geometry.euclidean.threed.Vector3D;
21  import org.hipparchus.util.FastMath;
22  
23  import java.lang.reflect.Field;
24  
25  import static org.junit.jupiter.api.Assertions.assertEquals;
26  import java.lang.reflect.InvocationTargetException;
27  import java.lang.reflect.Method;
28  
29  import org.junit.jupiter.api.Assertions;
30  import org.junit.jupiter.api.Test;
31  import org.orekit.bodies.GeodeticPoint;
32  import org.orekit.errors.OrekitException;
33  import org.orekit.rugged.api.AlgorithmId;
34  import org.orekit.rugged.errors.RuggedException;
35  import org.orekit.rugged.errors.RuggedMessages;
36  import org.orekit.rugged.intersection.AbstractAlgorithmTest;
37  import org.orekit.rugged.intersection.IntersectionAlgorithm;
38  import org.orekit.rugged.raster.CheckedPatternElevationUpdater;
39  import org.orekit.rugged.raster.Tile;
40  import org.orekit.rugged.raster.TileUpdater;
41  import org.orekit.rugged.raster.UpdatableTile;
42  import org.orekit.rugged.utils.ExtendedEllipsoid;
43  import org.orekit.rugged.utils.NormalizedGeodeticPoint;
44  
45  public class DuvenhageAlgorithmTest extends AbstractAlgorithmTest {
46  
47      protected IntersectionAlgorithm createAlgorithm(final TileUpdater updater, final int maxCachedTiles, final boolean isOverlappingTiles) {
48          return new DuvenhageAlgorithm(updater, maxCachedTiles, false, isOverlappingTiles);
49      }
50  
51      @Test
52      public void testNumericalIssueAtTileExit() {
53          setUpMayonVolcanoContext();
54          final IntersectionAlgorithm algorithm = createAlgorithm(updater, 8, true);
55          Vector3D position = new Vector3D(-3787079.6453602533, 5856784.405679551, 1655869.0582939098);
56          Vector3D los = new Vector3D( 0.5127552821932051, -0.8254313129088879, -0.2361041470463311);
57          GeodeticPoint intersection = algorithm.refineIntersection(earth, position, los,
58                                                                    algorithm.intersection(earth, position, los));
59          checkIntersection(position, los, intersection);
60      }
61  
62      @Test
63      public void testCrossingBeforeLineSegmentStart() {
64          setUpMayonVolcanoContext();
65          final IntersectionAlgorithm algorithm = createAlgorithm(updater, 8, true);
66          Vector3D position = new Vector3D(-3787079.6453602533, 5856784.405679551, 1655869.0582939098);
67          Vector3D los = new Vector3D( 0.42804005978915904, -0.8670291034054828, -0.2550338037664377);
68          GeodeticPoint intersection = algorithm.refineIntersection(earth, position, los,
69                                                                    algorithm.intersection(earth, position, los));
70          checkIntersection(position, los, intersection);
71      }
72  
73      @Test
74      public void testWrongPositionMissesGround() {
75          setUpMayonVolcanoContext();
76          final IntersectionAlgorithm algorithm = createAlgorithm(updater, 8, true);
77          Vector3D position = new Vector3D(7.551889113912788E9, -3.173692685491814E10, 1.5727517321541348E9);
78          Vector3D los = new Vector3D(0.010401349221417867, -0.17836068905951286, 0.9839101973923178);
79          try {
80              algorithm.intersection(earth, position, los);
81              Assertions.fail("an exception should have been thrown");
82          } catch (RuggedException re) {
83              Assertions.assertEquals(RuggedMessages.LINE_OF_SIGHT_DOES_NOT_REACH_GROUND, re.getSpecifier());
84          }
85      }
86  
87      @Test
88      public void testInconsistentTileUpdater() {
89          final int n = 1201;
90          final double size = FastMath.toRadians(1.0);
91          updater = new TileUpdater() {
92              public void updateTile(double latitude, double longitude, UpdatableTile tile) {
93                  
94                  double step = size / (n - 1);
95                  // this geometry is incorrect:
96                  // the specified latitude/longitude belong to rows/columns [1, n-1]
97                  // and not [0, n-2].
98                  tile.setGeometry(size * FastMath.floor(latitude / size) - 0.5 * step,
99                                   size * FastMath.floor(longitude / size) - 0.5 * step,
100                                  step, step, n, n);
101                 for (int i = 0; i < n; ++i) {
102                     for (int j = 0; j < n; ++j) {
103                         tile.setElevation(i, j, ((i + j) % 2 == 0) ? -7.0 : 224);
104                     }
105                 }
106             }
107         };
108         final IntersectionAlgorithm algorithm = createAlgorithm(updater, 8, true);
109         try {
110             algorithm.intersection(earth,
111                                    new Vector3D(-3010311.9672771087, 5307094.8081077365, 1852867.7919871407),
112                                    new Vector3D(0.3953329359154183, -0.8654901360032332, -0.30763402650162286));
113         } catch (RuggedException re) {
114             Assertions.assertEquals(RuggedMessages.TILE_WITHOUT_REQUIRED_NEIGHBORS_SELECTED, re.getSpecifier());
115         }
116     }
117 
118     @Test
119     public void testPureEastWestLOS() {
120         updater = new CheckedPatternElevationUpdater(FastMath.toRadians(1.0),1201, 41.0, 1563.0);
121         final IntersectionAlgorithm algorithm = createAlgorithm(updater, 8, true);
122         NormalizedGeodeticPoint gp =
123             algorithm.intersection(earth,
124                                    new Vector3D(-3041185.154503948, 6486750.132281409, -32335.022880173332),
125                                    new Vector3D(0.5660218606298548 , -0.8233939240951769, 0.040517885584811814));
126         Assertions.assertEquals(1164.35, gp.getAltitude(), 0.02);
127     }
128 
129     @Test
130     public void testParallelLOS() {
131         double size       = 0.125;
132         int    n          = 129;
133         double elevation1 = 0.0;
134         double elevation2 = 100.0;
135         updater = new CheckedPatternElevationUpdater(size, n, elevation1, elevation2);
136         MinMaxTreeTile northTile = new MinMaxTreeTileFactory().createTile();
137         updater.updateTile((3 * size) / 2, (3 * size) / 2, northTile);
138         MinMaxTreeTile southTile = new MinMaxTreeTileFactory().createTile();
139         updater.updateTile((-3 * size) / 2, (3 * size) / 2, southTile);
140         IntersectionAlgorithm algorithm = createAlgorithm(updater, 8, true);
141 
142         // line of sight in the South West corner
143         Assertions.assertEquals(northTile.getMinimumLongitude() - 0.0625 * northTile.getLongitudeStep(),
144                             findExit(algorithm, northTile,
145                                      new Vector3D( 6278799.86896170100,   788574.17965500170,   796074.70414069280),
146                                      new Vector3D( 0.09416282233912959,  0.01183204230132312, -0.99548649697728680)).getLongitude(),
147                                      1.0e-6);
148 
149         // line of sight in the West column
150         Assertions.assertEquals(northTile.getMinimumLongitude() - 0.0625 * northTile.getLongitudeStep(),
151                             findExit(algorithm, northTile,
152                                      new Vector3D(6278799.868961701, 788574.17965500171, 796074.7041406928),
153                                      new Vector3D(0.09231669916268806, 0.011600067441452849, -0.9956621241621375)).getLongitude(),
154                                      1.0e-6);
155 
156         // line of sight in the North-West corner
157         Assertions.assertEquals(northTile.getMinimumLongitude() - 0.0625 * northTile.getLongitudeStep(),
158                             findExit(algorithm, northTile,
159                                      new Vector3D( 6133039.79342824500,   770267.71434489540,  1568158.38266382620),
160                                      new Vector3D(-0.52028845147300570, -0.06537691642830394, -0.85148446025875800)).getLongitude(),
161                                      1.0e-6);
162         // line of sight in the North-East corner
163         Assertions.assertEquals(northTile.getMaximumLongitude() + 0.0625 * northTile.getLongitudeStep(),
164                             findExit(algorithm, northTile,
165                                      new Vector3D( 5988968.17708294100,  1529624.01701343130,  1568158.38266382620),
166                                      new Vector3D(-0.93877408645552440, -0.23970837882807683, -0.24747344851359457)).getLongitude(),
167                                      1.0e-6);
168 
169         // line of sight in the East column
170         Assertions.assertEquals(northTile.getMaximumLongitude() + 0.0625 * northTile.getLongitudeStep(),
171                             findExit(algorithm, northTile,
172                                      new Vector3D( 6106093.15406747100,  1559538.54861392200,   979886.66862965740),
173                                      new Vector3D(-0.18115090486319424, -0.04625542007869719,  0.98236693031707310)).getLongitude(),
174                                      1.0e-6);
175 
176         // line of sight in the South-East corner
177         Assertions.assertEquals(northTile.getMaximumLongitude() + 0.0625 * northTile.getLongitudeStep(),
178                             findExit(algorithm, northTile,
179                                      new Vector3D( 6131304.19368509600,  1565977.62301751650,   796074.70414069280),
180                                      new Vector3D( 0.09195297594530785,  0.02347944953986664, -0.99548649697728530)).getLongitude(),
181                                      1.0e-6);
182 
183         // line of sight in the South row
184         Assertions.assertEquals(northTile.getMinimumLatitude() - 0.0625 * northTile.getLatitudeStep(),
185                             findExit(algorithm, northTile,
186                                      new Vector3D(6251729.731998736, 984354.4963428857, 789526.5774750853),
187                                      new Vector3D(-0.15561499277355603, 0.9878177838164719, 0.0)).getLatitude(),
188                                      1.0e-6);
189 
190         // line of sight in the North row
191         Assertions.assertEquals(southTile.getMaximumLatitude() + 0.0625 * southTile.getLatitudeStep(),
192                             findExit(algorithm, southTile,
193                                      new Vector3D(6251729.731998736, 984354.4963428857, -789526.5774750853),
194                                      new Vector3D(-0.15561499277355603, 0.9878177838164719, 0.0)).getLatitude(),
195                                      1.0e-6);
196 
197     }
198     
199     @Test
200     public void testAlgorithmId() {
201         setUpMayonVolcanoContext();
202  
203         final IntersectionAlgorithm algorithm = new DuvenhageAlgorithm(updater, 8, false, true);
204         assertEquals(AlgorithmId.DUVENHAGE, algorithm.getAlgorithmId());
205         
206         final IntersectionAlgorithm algorithmFlatBody = new DuvenhageAlgorithm(updater, 8, true, true);
207         assertEquals(AlgorithmId.DUVENHAGE_FLAT_BODY, algorithmFlatBody.getAlgorithmId());
208     }
209 
210     private NormalizedGeodeticPoint findExit(IntersectionAlgorithm algorithm, Tile tile, Vector3D position, Vector3D los) {
211 
212         try {
213             Method findExit = DuvenhageAlgorithm.class.getDeclaredMethod("findExit",
214                                                                          Tile.class,
215                                                                          ExtendedEllipsoid.class,
216                                                                          Vector3D.class, Vector3D.class);
217             findExit.setAccessible(true);
218             Object limitPoint = findExit.invoke(algorithm, tile, earth, position, los);
219             Class<?> limitPointCls = null;
220             for (Class<?> c : DuvenhageAlgorithm.class.getDeclaredClasses()) {
221                 if (c.getName().endsWith("LimitPoint")) {
222                     limitPointCls = c;
223                 }
224             }
225             Field pointField = limitPointCls.getDeclaredField("point");
226             pointField.setAccessible(true);
227             return (NormalizedGeodeticPoint) pointField.get(limitPoint);
228 
229         } catch (NoSuchMethodException nsme) {
230             Assertions.fail(nsme.getLocalizedMessage());
231         } catch (SecurityException se) {
232             Assertions.fail(se.getLocalizedMessage());
233         } catch (IllegalAccessException iae) {
234             Assertions.fail(iae.getLocalizedMessage());
235         } catch (IllegalArgumentException iae) {
236             Assertions.fail(iae.getLocalizedMessage());
237         } catch (NoSuchFieldException nsfe) {
238             Assertions.fail(nsfe.getLocalizedMessage());
239         } catch (InvocationTargetException e) {
240             if (e.getCause() instanceof RuggedException) {
241                 throw (RuggedException) e.getCause();
242             } else {
243                 throw (OrekitException) e.getCause();
244             }
245         }
246         return null;
247     }
248 
249 }