1   /* Copyright 2002-2022 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.geometry.fov;
18  
19  import java.util.List;
20  
21  import org.hipparchus.geometry.euclidean.threed.RotationOrder;
22  import org.hipparchus.geometry.euclidean.threed.Vector3D;
23  import org.hipparchus.util.FastMath;
24  import org.junit.After;
25  import org.junit.Assert;
26  import org.junit.Before;
27  import org.junit.Test;
28  import org.orekit.Utils;
29  import org.orekit.attitudes.AttitudeProvider;
30  import org.orekit.attitudes.LofOffset;
31  import org.orekit.attitudes.NadirPointing;
32  import org.orekit.bodies.GeodeticPoint;
33  import org.orekit.bodies.OneAxisEllipsoid;
34  import org.orekit.errors.OrekitException;
35  import org.orekit.errors.OrekitMessages;
36  import org.orekit.frames.FramesFactory;
37  import org.orekit.frames.LOFType;
38  import org.orekit.frames.TopocentricFrame;
39  import org.orekit.frames.Transform;
40  import org.orekit.geometry.fov.PolygonalFieldOfView.DefiningConeType;
41  import org.orekit.orbits.KeplerianOrbit;
42  import org.orekit.orbits.Orbit;
43  import org.orekit.propagation.Propagator;
44  import org.orekit.propagation.SpacecraftState;
45  import org.orekit.propagation.analytical.KeplerianPropagator;
46  import org.orekit.propagation.events.VisibilityTrigger;
47  import org.orekit.time.AbsoluteDate;
48  import org.orekit.utils.Constants;
49  import org.orekit.utils.IERSConventions;
50  import org.orekit.utils.PVCoordinates;
51  
52  public class PolygonalFieldOfViewTest {
53  
54      @Test
55      public void testRegularPolygon() {
56          double delta          = 0.25;
57          double margin         = 0.01;
58          double maxAreaError   = 0;
59          double maxOffsetError = 0;
60          for (int n = 3; n < 32; ++n) {
61              PolygonalFieldOfView base = new PolygonalFieldOfView(Vector3D.PLUS_K,
62                                                                   DefiningConeType.INSIDE_CONE_TOUCHING_POLYGON_AT_EDGES_MIDDLE,
63                                                                   Vector3D.PLUS_I, delta, n, margin);
64              PolygonalFieldOfView fov  = new PolygonalFieldOfView(base.getZone(), margin);
65              double eta = FastMath.acos(FastMath.sin(FastMath.PI / n) * FastMath.cos(delta));
66              double theoreticalArea = 2 * n * eta - (n - 2) * FastMath.PI;
67              double areaError = theoreticalArea - fov.getZone().getSize();
68              maxAreaError = FastMath.max(FastMath.abs(areaError), maxAreaError);
69              for (double lambda = -0.5 * FastMath.PI; lambda < 0.5 * FastMath.PI; lambda += 0.1) {
70                  Vector3D v = new Vector3D(0.0, lambda).scalarMultiply(1.0e6);
71                  double theoreticalOffset = 0.5 * FastMath.PI - lambda - delta - margin;
72                  double offset = fov.offsetFromBoundary(v, 0.0, VisibilityTrigger.VISIBLE_ONLY_WHEN_FULLY_IN_FOV);
73                  if (theoreticalOffset > 0.01) {
74                      // the offsetFromBoundary method may use the fast approximate
75                      // method, so we cannot check the error accurately
76                      // we know however that the fast method will underestimate the offset
77  
78                      Assert.assertTrue(offset > 0);
79                      Assert.assertTrue(offset <= theoreticalOffset + 5e-16);
80                  } else {
81                      double offsetError = theoreticalOffset - offset;
82                      maxOffsetError = FastMath.max(FastMath.abs(offsetError), maxOffsetError);
83                  }
84                  Assert.assertEquals(-margin,
85                                      fov.offsetFromBoundary(fov.projectToBoundary(v), 0.0, VisibilityTrigger.VISIBLE_ONLY_WHEN_FULLY_IN_FOV),
86                                      1.0e-12);
87              }
88          }
89          Assert.assertEquals(0.0, maxAreaError,   5.0e-14);
90          Assert.assertEquals(0.0, maxOffsetError, 2.0e-15);
91      }
92  
93      @Test
94      public void testNoFootprintInside() {
95          Utils.setDataRoot("regular-data");
96          PolygonalFieldOfView fov = new PolygonalFieldOfView(Vector3D.PLUS_K,
97                                                              DefiningConeType.INSIDE_CONE_TOUCHING_POLYGON_AT_EDGES_MIDDLE,
98                                                              Vector3D.PLUS_I,
99                                                              FastMath.toRadians(3.0), 6, 0.0);
100         OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
101                                                       Constants.WGS84_EARTH_FLATTENING,
102                                                       FramesFactory.getITRF(IERSConventions.IERS_2010, true));
103         Transform fovToBody   = new Transform(AbsoluteDate.J2000_EPOCH, new Vector3D(5e6, 3e6, 2e6));
104         try {
105             fov.getFootprint(fovToBody, earth, FastMath.toRadians(0.1));
106             Assert.fail("an exception should have been thrown");
107         } catch (OrekitException oe) {
108             Assert.assertEquals(OrekitMessages.POINT_INSIDE_ELLIPSOID, oe.getSpecifier());
109         }
110     }
111 
112     @Test
113     public void testNadirHexagonalFootprint() {
114         doTest(new PolygonalFieldOfView(Vector3D.PLUS_K,
115                                         DefiningConeType.INSIDE_CONE_TOUCHING_POLYGON_AT_EDGES_MIDDLE,
116                                         Vector3D.PLUS_I,
117                                         FastMath.toRadians(3.0), 6, 0.0),
118                new NadirPointing(orbit.getFrame(), earth),
119                210, 84.6497, 85.3729, 181052.2, 209092.8);
120     }
121 
122     @Test
123     public void testRollPitchYawHexagonalFootprint() {
124         doTest(new PolygonalFieldOfView(Vector3D.PLUS_K,
125                                         DefiningConeType.INSIDE_CONE_TOUCHING_POLYGON_AT_EDGES_MIDDLE,
126                                         Vector3D.PLUS_I,
127                                         FastMath.toRadians(3.0), 6, 0.0),
128                new LofOffset(orbit.getFrame(), LOFType.LVLH_CCSDS, RotationOrder.XYZ,
129                              FastMath.toRadians(10),
130                              FastMath.toRadians(20),
131                              FastMath.toRadians(5)),
132                210, 48.0026, 60.1975, 1221543.6, 1804921.6);
133     }
134 
135     @Test
136     public void testFOVPartiallyTruncatedAtLimb() {
137         doTest(new PolygonalFieldOfView(Vector3D.PLUS_K,
138                                         DefiningConeType.INSIDE_CONE_TOUCHING_POLYGON_AT_EDGES_MIDDLE,
139                                         Vector3D.PLUS_I,
140                                         FastMath.toRadians(40.0), 6, 0.0),
141                new NadirPointing(orbit.getFrame(), earth),
142                2448, 0.0, 7.9089, 4583054.6, 5347029.8);
143     }
144 
145     @Test
146     public void testFOVLargerThanEarth() {
147         doTest(new PolygonalFieldOfView(Vector3D.PLUS_K,
148                                         DefiningConeType.INSIDE_CONE_TOUCHING_POLYGON_AT_EDGES_MIDDLE,
149                                         Vector3D.PLUS_I,
150                                         FastMath.toRadians(45.0), 6, 0.0),
151                new NadirPointing(orbit.getFrame(), earth),
152                2337, 0.0, 0.0, 5323032.8, 5347029.8);
153     }
154 
155     @Test
156     public void testFOVLargerThanEarthOld() {
157         Utils.setDataRoot("regular-data");
158         PolygonalFieldOfView fov = new PolygonalFieldOfView(Vector3D.PLUS_K,
159                                                             DefiningConeType.INSIDE_CONE_TOUCHING_POLYGON_AT_EDGES_MIDDLE,
160                                                             Vector3D.PLUS_I,
161                                                             FastMath.toRadians(45.0), 6, 0.0);
162         OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
163                                                       Constants.WGS84_EARTH_FLATTENING,
164                                                       FramesFactory.getITRF(IERSConventions.IERS_2010, true));
165         KeplerianOrbit orbit = new KeplerianOrbit(new PVCoordinates(new Vector3D(7.0e6, 1.0e6, 4.0e6),
166                                                                     new Vector3D(-500.0, 8000.0, 1000.0)),
167                                                   FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH,
168                                                   Constants.EIGEN5C_EARTH_MU);
169         Propagator propagator = new KeplerianPropagator(orbit);
170         propagator.setAttitudeProvider(new NadirPointing(orbit.getFrame(), earth));
171         SpacecraftState state = propagator.propagate(orbit.getDate().shiftedBy(1000.0));
172         Transform inertToBody = state.getFrame().getTransformTo(earth.getBodyFrame(), state.getDate());
173         Transform fovToBody   = new Transform(state.getDate(),
174                                               state.toTransform().getInverse(),
175                                               inertToBody);
176         List<List<GeodeticPoint>> footprint = fov.getFootprint(fovToBody, earth, FastMath.toRadians(1.0));
177         Vector3D subSat = earth.projectToGround(state.getPVCoordinates(earth.getBodyFrame()).getPosition(),
178                                                 state.getDate(), earth.getBodyFrame());
179         Assert.assertEquals(1, footprint.size());
180         List<GeodeticPoint> loop = footprint.get(0);
181         Assert.assertEquals(234, loop.size());
182         double minEl   = Double.POSITIVE_INFINITY;
183         double maxEl = 0;
184         double minDist = Double.POSITIVE_INFINITY;
185         double maxDist = 0;
186         for (int i = 0; i < loop.size(); ++i) {
187             Assert.assertEquals(0.0, loop.get(i).getAltitude(), 3.0e-7);
188             TopocentricFrame topo = new TopocentricFrame(earth, loop.get(i), "atLimb");
189             final double elevation = topo.getElevation(state.getPVCoordinates().getPosition(),
190                                                        state.getFrame(), state.getDate());
191             minEl = FastMath.min(minEl, elevation);
192             maxEl = FastMath.max(maxEl, elevation);
193             final double dist = Vector3D.distance(subSat, earth.transform(loop.get(i)));
194             minDist = FastMath.min(minDist, dist);
195             maxDist = FastMath.max(maxDist, dist);
196         }
197         Assert.assertEquals(0.0,       FastMath.toDegrees(minEl), 2.0e-12);
198         Assert.assertEquals(0.0,       FastMath.toDegrees(maxEl), 1.7e-12);
199         Assert.assertEquals(5323036.6, minDist, 1.0);
200         Assert.assertEquals(5347029.8, maxDist, 1.0);
201     }
202 
203     @Test
204     public void testFOVAwayFromEarth() {
205         PolygonalFieldOfView fov = new PolygonalFieldOfView(Vector3D.MINUS_K,
206                                                             DefiningConeType.INSIDE_CONE_TOUCHING_POLYGON_AT_EDGES_MIDDLE,
207                                                             Vector3D.PLUS_I,
208                                                             FastMath.toRadians(3.0), 6, 0.0);
209         Propagator propagator = new KeplerianPropagator(orbit);
210         propagator.setAttitudeProvider(new NadirPointing(orbit.getFrame(), earth));
211         SpacecraftState state = propagator.propagate(orbit.getDate().shiftedBy(1000.0));
212         Transform inertToBody = state.getFrame().getTransformTo(earth.getBodyFrame(), state.getDate());
213         Transform fovToBody   = new Transform(state.getDate(),
214                                               state.toTransform().getInverse(),
215                                               inertToBody);
216         List<List<GeodeticPoint>> footprint = fov.getFootprint(fovToBody, earth, FastMath.toRadians(1.0));
217         Assert.assertEquals(0, footprint.size());
218     }
219 
220     private void doTest(final PolygonalFieldOfView fov, final AttitudeProvider attitude, final int expectedPoints,
221                         final double expectedMinElevation, final double expectedMaxElevation,
222                         final double expectedMinDist, final double expectedMaxDist) {
223 
224         Propagator propagator = new KeplerianPropagator(orbit);
225         propagator.setAttitudeProvider(attitude);
226         SpacecraftState state = propagator.propagate(orbit.getDate().shiftedBy(1000.0));
227         Transform inertToBody = state.getFrame().getTransformTo(earth.getBodyFrame(), state.getDate());
228         Transform fovToBody   = new Transform(state.getDate(),
229                                               state.toTransform().getInverse(),
230                                               inertToBody);
231         List<List<GeodeticPoint>> footprint = fov.getFootprint(fovToBody, earth, FastMath.toRadians(0.1));
232         Vector3D subSat = earth.projectToGround(state.getPVCoordinates(earth.getBodyFrame()).getPosition(),
233                                                 state.getDate(), earth.getBodyFrame());
234         Assert.assertEquals(1, footprint.size());
235         List<GeodeticPoint> loop = footprint.get(0);
236         Assert.assertEquals(expectedPoints, loop.size());
237         double minEl     = Double.POSITIVE_INFINITY;
238         double maxEl     = 0;
239         double minDist   = Double.POSITIVE_INFINITY;
240         double maxDist   = 0;
241         for (int i = 0; i < loop.size(); ++i) {
242 
243             Assert.assertEquals(0.0, loop.get(i).getAltitude(), 9.0e-9);
244 
245             TopocentricFrame topo = new TopocentricFrame(earth, loop.get(i), "onFootprint");
246             final double elevation = topo.getElevation(state.getPVCoordinates().getPosition(),
247                                                        state.getFrame(), state.getDate());
248             if (elevation > 0.001) {
249                 Vector3D los = fovToBody.getInverse().transformPosition(earth.transform(loop.get(i)));
250                 Assert.assertEquals(-fov.getMargin(),
251                                     fov.offsetFromBoundary(los, 0.0, VisibilityTrigger.VISIBLE_ONLY_WHEN_FULLY_IN_FOV),
252                                     4.0e-15);
253                 Assert.assertEquals(0.125 - fov.getMargin(),
254                                     fov.offsetFromBoundary(los, 0.125, VisibilityTrigger.VISIBLE_ONLY_WHEN_FULLY_IN_FOV),
255                                     4.0e-15);
256                 Assert.assertEquals(-0.125 - fov.getMargin(),
257                                     fov.offsetFromBoundary(los, 0.125, VisibilityTrigger.VISIBLE_AS_SOON_AS_PARTIALLY_IN_FOV),
258                                     4.0e-15);
259             }
260             minEl = FastMath.min(minEl, elevation);
261             maxEl = FastMath.max(maxEl, elevation);
262             final double dist = Vector3D.distance(subSat, earth.transform(loop.get(i)));
263             minDist = FastMath.min(minDist, dist);
264             maxDist = FastMath.max(maxDist, dist);
265 
266         }
267 
268         Assert.assertEquals(expectedMinElevation, FastMath.toDegrees(minEl), 0.001);
269         Assert.assertEquals(expectedMaxElevation, FastMath.toDegrees(maxEl), 0.001);
270         Assert.assertEquals(expectedMinDist,      minDist,                   1.0);
271         Assert.assertEquals(expectedMaxDist,      maxDist,                   1.0);
272 
273     }
274 
275     @Before
276     public void setUp() {
277         Utils.setDataRoot("regular-data");
278         earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
279                                      Constants.WGS84_EARTH_FLATTENING,
280                                      FramesFactory.getITRF(IERSConventions.IERS_2010, true));
281         orbit = new KeplerianOrbit(new PVCoordinates(new Vector3D(7.0e6, 1.0e6, 4.0e6),
282                                                      new Vector3D(-500.0, 8000.0, 1000.0)),
283                                    FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH,
284                                    Constants.EIGEN5C_EARTH_MU);
285     }
286 
287     @After
288     public void tearDown() {
289         earth = null;
290         orbit = null;
291     }
292 
293     private OneAxisEllipsoid earth;
294     private Orbit            orbit;
295     
296 }