FootprintOverlapDetector.java

  1. /* Copyright 2002-2019 CS Systèmes d'Information
  2.  * Licensed to CS Systèmes d'Information (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.propagation.events;

  18. import java.util.ArrayList;
  19. import java.util.List;

  20. import org.hipparchus.geometry.enclosing.EnclosingBall;
  21. import org.hipparchus.geometry.euclidean.threed.Vector3D;
  22. import org.hipparchus.geometry.spherical.twod.Edge;
  23. import org.hipparchus.geometry.spherical.twod.S2Point;
  24. import org.hipparchus.geometry.spherical.twod.Sphere2D;
  25. import org.hipparchus.geometry.spherical.twod.SphericalPolygonsSet;
  26. import org.hipparchus.geometry.spherical.twod.Vertex;
  27. import org.hipparchus.ode.events.Action;
  28. import org.hipparchus.util.FastMath;
  29. import org.orekit.bodies.BodyShape;
  30. import org.orekit.bodies.GeodeticPoint;
  31. import org.orekit.bodies.OneAxisEllipsoid;
  32. import org.orekit.frames.Transform;
  33. import org.orekit.models.earth.tessellation.DivertedSingularityAiming;
  34. import org.orekit.models.earth.tessellation.EllipsoidTessellator;
  35. import org.orekit.propagation.SpacecraftState;
  36. import org.orekit.propagation.events.handlers.EventHandler;
  37. import org.orekit.propagation.events.handlers.StopOnIncreasing;

  38. /** Detector triggered by geographical region entering/leaving a spacecraft sensor
  39.  * {@link FieldOfView Field Of View}.
  40.  * <p>
  41.  * This detector is a mix between to {@link FieldOfViewDetector} and {@link
  42.  * GeographicZoneDetector}. Similar to the first detector above, it triggers events
  43.  * related to entry/exit of targets in a Field Of View, taking attitude into account.
  44.  * Similar to the second detector above, its target is an entire geographic region
  45.  * (which can even be split in several non-connected patches and can have holes).
  46.  * </p>
  47.  * <p>
  48.  * This detector is typically used for ground observation missions with agile
  49.  * satellites than can look away from nadir.
  50.  * </p>
  51.  * <p>The default implementation behavior is to {@link Action#CONTINUE continue}
  52.  * propagation at FOV entry and to {@link Action#STOP stop} propagation
  53.  * at FOV exit. This can be changed by calling
  54.  * {@link #withHandler(EventHandler)} after construction.</p>
  55.  * @see org.orekit.propagation.Propagator#addEventDetector(EventDetector)
  56.  * @see FieldOfViewDetector
  57.  * @see GeographicZoneDetector
  58.  * @author Luc Maisonobe
  59.  * @since 7.1
  60.  */
  61. public class FootprintOverlapDetector extends AbstractDetector<FootprintOverlapDetector> {

  62.     /** Field of view. */
  63.     private final FieldOfView fov;

  64.     /** Body on which the geographic zone is defined. */
  65.     private final OneAxisEllipsoid body;

  66.     /** Geographic zone to consider. */
  67.     private final SphericalPolygonsSet zone;

  68.     /** Linear step used for sampling the geographic zone. */
  69.     private final double samplingStep;

  70.     /** Sampling of the geographic zone. */
  71.     private final List<SamplingPoint> sampledZone;

  72.     /** Center of the spherical cap surrounding the zone. */
  73.     private final Vector3D capCenter;

  74.     /** Cosine of the radius of the spherical cap surrounding the zone. */
  75.     private final double capCos;

  76.     /** Sine of the radius of the spherical cap surrounding the zone. */
  77.     private final double capSin;

  78.     /** Build a new instance.
  79.      * <p>The maximal interval between distance to FOV boundary checks should
  80.      * be smaller than the half duration of the minimal pass to handle,
  81.      * otherwise some short passes could be missed.</p>
  82.      * @param fov sensor field of view
  83.      * @param body body on which the geographic zone is defined
  84.      * @param zone geographic zone to consider
  85.      * @param samplingStep linear step used for sampling the geographic zone (in meters)
  86.      */
  87.     public FootprintOverlapDetector(final FieldOfView fov,
  88.                                     final OneAxisEllipsoid body,
  89.                                     final SphericalPolygonsSet zone,
  90.                                     final double samplingStep) {
  91.         this(DEFAULT_MAXCHECK, DEFAULT_THRESHOLD, DEFAULT_MAX_ITER,
  92.              new StopOnIncreasing<FootprintOverlapDetector>(),
  93.              fov, body, zone, samplingStep, sample(body, zone, samplingStep));
  94.     }

  95.     /** Private constructor with full parameters.
  96.      * <p>
  97.      * This constructor is private as users are expected to use the builder
  98.      * API with the various {@code withXxx()} methods to set up the instance
  99.      * in a readable manner without using a huge amount of parameters.
  100.      * </p>
  101.      * @param maxCheck maximum checking interval (s)
  102.      * @param threshold convergence threshold (s)
  103.      * @param maxIter maximum number of iterations in the event time search
  104.      * @param handler event handler to call at event occurrences
  105.      * @param body body on which the geographic zone is defined
  106.      * @param zone geographic zone to consider
  107.      * @param fov sensor field of view
  108.      * @param sampledZone sampling of the geographic zone
  109.      * @param samplingStep linear step used for sampling the geographic zone (in meters)
  110.      */
  111.     private FootprintOverlapDetector(final double maxCheck, final double threshold,
  112.                                      final int maxIter, final EventHandler<? super FootprintOverlapDetector> handler,
  113.                                      final FieldOfView fov,
  114.                                      final OneAxisEllipsoid body,
  115.                                      final SphericalPolygonsSet zone,
  116.                                      final double samplingStep,
  117.                                      final List<SamplingPoint> sampledZone) {

  118.         super(maxCheck, threshold, maxIter, handler);
  119.         this.fov          = fov;
  120.         this.body         = body;
  121.         this.samplingStep = samplingStep;
  122.         this.zone         = zone;
  123.         this.sampledZone  = sampledZone;

  124.         final EnclosingBall<Sphere2D, S2Point> cap = zone.getEnclosingCap();
  125.         this.capCenter    = cap.getCenter().getVector();
  126.         this.capCos       = FastMath.cos(cap.getRadius());
  127.         this.capSin       = FastMath.sin(cap.getRadius());

  128.     }

  129.     /** Sample the region.
  130.      * @param body body on which the geographic zone is defined
  131.      * @param zone geographic zone to consider
  132.      * @param samplingStep  linear step used for sampling the geographic zone (in meters)
  133.      * @return sampling points
  134.      */
  135.     private static List<SamplingPoint> sample(final OneAxisEllipsoid body,
  136.                                               final SphericalPolygonsSet zone,
  137.                                               final double samplingStep) {

  138.         final List<SamplingPoint> sampledZone = new ArrayList<SamplingPoint>();

  139.         // sample the zone boundary
  140.         final List<Vertex> boundary = zone.getBoundaryLoops();
  141.         for (final Vertex loopStart : boundary) {
  142.             int count = 0;
  143.             for (Vertex v = loopStart; count == 0 || v != loopStart; v = v.getOutgoing().getEnd()) {
  144.                 ++count;
  145.                 final Edge edge = v.getOutgoing();
  146.                 final int n = (int) FastMath.ceil(edge.getLength() * body.getEquatorialRadius() / samplingStep);
  147.                 for (int i = 0; i < n; ++i) {
  148.                     final S2Point intermediate = new S2Point(edge.getPointAt(i * edge.getLength() / n));
  149.                     final GeodeticPoint gp = new GeodeticPoint(0.5 * FastMath.PI - intermediate.getPhi(),
  150.                                                                intermediate.getTheta(), 0.0);
  151.                     sampledZone.add(new SamplingPoint(body.transform(gp), gp.getZenith()));
  152.                 }
  153.             }
  154.         }

  155.         // sample the zone interior
  156.         final EllipsoidTessellator tessellator =
  157.                         new EllipsoidTessellator(body, new DivertedSingularityAiming(zone), 1);
  158.         final List<List<GeodeticPoint>> gpSample = tessellator.sample(zone, samplingStep, samplingStep);
  159.         for (final List<GeodeticPoint> list : gpSample) {
  160.             for (final GeodeticPoint gp : list) {
  161.                 sampledZone.add(new SamplingPoint(body.transform(gp), gp.getZenith()));
  162.             }
  163.         }

  164.         return sampledZone;

  165.     }

  166.     /** {@inheritDoc} */
  167.     @Override
  168.     protected FootprintOverlapDetector create(final double newMaxCheck, final double newThreshold,
  169.                                               final int newMaxIter,
  170.                                               final EventHandler<? super FootprintOverlapDetector> newHandler) {
  171.         return new FootprintOverlapDetector(newMaxCheck, newThreshold, newMaxIter, newHandler,
  172.                                             fov, body, zone, samplingStep, sampledZone);
  173.     }

  174.     /** Get the geographic zone triggering the events.
  175.      * <p>
  176.      * The zone is mapped on the unit sphere
  177.      * </p>
  178.      * @return geographic zone triggering the events
  179.      */
  180.     public SphericalPolygonsSet getZone() {
  181.         return zone;
  182.     }

  183.     /** Get the Field Of View.
  184.      * @return Field Of View
  185.      */
  186.     public FieldOfView getFieldOfView() {
  187.         return fov;
  188.     }

  189.     /** Get the body on which the geographic zone is defined.
  190.      * @return body on which the geographic zone is defined
  191.      */
  192.     public BodyShape getBody() {
  193.         return body;
  194.     }

  195.     /** {@inheritDoc}
  196.      * <p>
  197.      * The g function value is the minimum offset among the region points
  198.      * with respect to the Field Of View boundary. It is positive if all region
  199.      * points are outside of the Field Of View, and negative if at least some
  200.      * of the region points are inside of the Field Of View. The minimum is
  201.      * computed by sampling the region, considering only the points for which
  202.      * the spacecraft is above the horizon. The accuracy of the detection
  203.      * depends on the linear sampling step set at detector construction. If
  204.      * the spacecraft is below horizon for all region points, an arbitrary
  205.      * positive value is returned.
  206.      * </p>
  207.      * <p>
  208.      * As per the previous definition, when the region enters the Field Of
  209.      * View, a decreasing event is generated, and when the region leaves
  210.      * the Field Of View, an increasing event is generated.
  211.      * </p>
  212.      */
  213.     public double g(final SpacecraftState s) {

  214.         // initial arbitrary positive value
  215.         double value = FastMath.PI;

  216.         // get spacecraft position in body frame
  217.         final Vector3D      scBody      = s.getPVCoordinates(body.getBodyFrame()).getPosition();

  218.         // map the point to a sphere
  219.         final GeodeticPoint gp          = body.transform(scBody, body.getBodyFrame(), s.getDate());
  220.         final S2Point       s2p         = new S2Point(gp.getLongitude(), 0.5 * FastMath.PI - gp.getLatitude());

  221.         // for faster computation, we start using only the surrounding cap, to filter out
  222.         // far away points (which correspond to most of the points if the zone is small)
  223.         final Vector3D p   = s2p.getVector();
  224.         final double   dot = Vector3D.dotProduct(p, capCenter);
  225.         if (dot < capCos) {
  226.             // the spacecraft is outside of the cap, look for the closest cap point
  227.             final Vector3D t     = p.subtract(dot, capCenter).normalize();
  228.             final Vector3D close = new Vector3D(capCos, capCenter, capSin, t);
  229.             if (Vector3D.dotProduct(p, close) < -0.01) {
  230.                 // the spacecraft is not visible from the cap edge,
  231.                 // even taking some margin into account for sphere/ellipsoid different shapes
  232.                 // this induces no points in the zone can see the spacecraft,
  233.                 // we can return the arbitrary initial positive value without performing further computation
  234.                 return value;
  235.             }
  236.         }

  237.         // the spacecraft may be visible from some points in the zone, check them all
  238.         final Transform bodyToSc = new Transform(s.getDate(),
  239.                                                  body.getBodyFrame().getTransformTo(s.getFrame(), s.getDate()),
  240.                                                  s.toTransform());
  241.         for (final SamplingPoint point : sampledZone) {
  242.             final Vector3D lineOfSightBody = point.getPosition().subtract(scBody);
  243.             if (Vector3D.dotProduct(lineOfSightBody, point.getZenith()) <= 0) {
  244.                 // spacecraft is above this sample point local horizon
  245.                 // get line of sight in spacecraft frame
  246.                 final double offset = fov.offsetFromBoundary(bodyToSc.transformVector(lineOfSightBody));
  247.                 value = FastMath.min(value, offset);
  248.             }
  249.         }

  250.         return value;

  251.     }

  252.     /** Container for sampling points. */
  253.     private static class SamplingPoint {

  254.         /** Position of the point. */
  255.         private final Vector3D position;

  256.         /** Zenith vector of the point. */
  257.         private final Vector3D zenith;

  258.         /** Simple constructor.
  259.          * @param position position of the point
  260.          * @param zenith zenith vector of the point
  261.          */
  262.         SamplingPoint(final Vector3D position, final Vector3D zenith) {
  263.             this.position = position;
  264.             this.zenith   = zenith;
  265.         }

  266.         /** Get the point position.
  267.          * @return point position
  268.          */
  269.         public Vector3D getPosition() {
  270.             return position;
  271.         }

  272.         /** Get the point zenith vector.
  273.          * @return point zenith vector
  274.          */
  275.         public Vector3D getZenith() {
  276.             return zenith;
  277.         }

  278.     }

  279. }