Define a ground area target

1

Hello, I have encountered an issue. I would like to define a ground area target, such as a rectangular ground area. Could you please explain how to use the Orekit library for this?

2

Hi @YangYang

You can try using the EllipsoidTessellator.buildSimpleZone() method. This method requires a list of GeodeticPoint objects defining the ground area. Please find below an example with 4 points:

GeodeticPoint point1 = new GeodeticPoint(lat[0],lon[0],alt[0])
GeodeticPoint point2 = new GeodeticPoint(lat[1],lon[1],alt[1])
GeodeticPoint point3 = new GeodeticPoint(lat[2],lon[2],alt[2])
GeodeticPoint point4 = new GeodeticPoint(lat[3],lon[3],alt[3])
SphericalPolygonsSet region = EllipsoidTessellator.buildSimpleZone(1.0e-10, point1,point2,point3,point4)

(1) it is important to note that lat and lon shall be in radians. alt shall be in meters.
(2) This method returns a SphericalPolygonsSet that you can use for diverse applications. Like inside a GeographicZoneDetector if you want to detect when a satellite is visible from the area.

Best regards,
Bryan

1 Like
3

Thank you for your answer. Wishing you a pleasant work and life！

4

So, if I want to define a pentagon, hexagon, heptagon… what exactly should I do?

5

You can use 5, 6, or 7 GeodeticPoint objects.

6

Thank you very much！

7

Another comment: take care of the order of the points. The interior of the region is considered to be on the left hand side when you travel from one point to the next. So two regions defined by the exact same list of points but in reversed order are complementary to each other.

9

How do we understand this complementarity? Is it like the left handed helix rule and the right handed helix rule in physics?

10

Hello, I have defined a ground area target and used the GeographicZoneDetector to define a detector. Now I want to use the detector to observe the ground target. Which method in the orekit library should I call? Or what is the specific usage code?

11

Not exactly, but close.
As the sphere is a finite closed manifold, when we define a boundary separating two regions, there is no canonical way to say one is the interior and one is the exterior as we could do in an infinite 2D plane where interior would be the finite area and exterior the infinite area. So on the 2D sphere, we have to use come conventions. The one we adopted is that the boundary is oriented and if someone travels along the boundary according to this orientation, the interior is on the left hand side and the exterior is on the right hand side. Complementary here means that the union of the two regions is the whole sphere.

12

This is a different problem. GeographicZoneDetector is used to define when satellite hovers above a complete zone, not when the satellite can observe a single point target.

If what you want is observing a point from a satellite field of view, the detector to use is FieldOfViewDetector and then the zone represents the field of view from the satellite perspective.

If you need both a zone on ground and a field of view on board, i.e. when part of the satellite field of view can see the zone, even if the satellite is not directly above it (for example because the field of view is wide enough or the satellite is pointing away from nadir) then you need to use FootprintOverlapDetector. Beware however that this detector is far from perfect. It relies on sampling the ground zone in a number of points, and checks all these points with respect to the on-board field of view. If the ground zone is large (or the sampling step is small), then it is really computationally costly to do that. This is a known limitation of this detector. We want to fix that, but did not find time to address this yet.

13

So, the direction that the left hand can wrap around (clockwise) is the inside, and the direction that the right hand can wrap around (counterclockwise) is the outside?

14

No, it is the other way round.

15

One way to check is to call size() on the spherical region. It will return the surface area in steradians. If you define a “small” zone (say an island) the area should be smaller that 2π steradians (i.e. less than one hemisphere). If the area is larger than that, then you have defined the complementary region, i.e. the whole Earth minus the island.

1 Like
16

I understand. This is like calculating the size of a trigonometric function. If we ignore the difference between obtuse and acute angles, sin (85deg)=sin (180deg-85deg).
Thank you very much for your help. Wishing you a pleasant work and life

1 Like
17

Hello, I have encountered some problems again now. I have defined a rectangular observation area on the ground, and I have defined a satellite for observation. But I found that the execution results of the observation were not the same as the observation results of STK. May I ask if the observation task performed here is the calculation of the visible time period of the satellite in the ground area or the time period of the satellite’s sub satellite point trajectory in the ground area?
Best wishes!
import java.io.File;
import java.util.Locale;
import org.hipparchus.geometry.euclidean.threed.Vector3D;
import org.hipparchus.geometry.spherical.twod.SphericalPolygonsSet;
import org.hipparchus.ode.events.Action;
import org.hipparchus.util.FastMath;
import org.orekit.bodies.BodyShape;
import org.orekit.bodies.GeodeticPoint;
import org.orekit.bodies.OneAxisEllipsoid;
import org.orekit.data.DataContext;
import org.orekit.data.DataProvidersManager;
import org.orekit.data.DirectoryCrawler;
import org.orekit.errors.OrekitException;
import org.orekit.frames.Frame;
import org.orekit.frames.FramesFactory;
import org.orekit.frames.TopocentricFrame;
import org.orekit.models.earth.tessellation.EllipsoidTessellator;
import org.orekit.orbits.KeplerianOrbit;
import org.orekit.orbits.Orbit;
import org.orekit.propagation.Propagator;
import org.orekit.propagation.SpacecraftState;
import org.orekit.propagation.analytical.KeplerianPropagator;
import org.orekit.propagation.events.ElevationDetector;
import org.orekit.propagation.events.EventDetector;
import org.orekit.propagation.events.GeographicZoneDetector;
import org.orekit.time.AbsoluteDate;
import org.orekit.time.TimeScalesFactory;
import org.orekit.utils.Constants;
import org.orekit.utils.IERSConventions;
import org.orekit.utils.PVCoordinates;

/** Orekit tutorial for special event detection.

  • This tutorial shows how to easily check for visibility between a satellite and a ground station.

  • @author Pascal Parraud
    */
    public class RectangularRegionOnEarth {

    /** Private constructor for utility class. */
    private RectangularRegionOnEarth() {
    // empty
    }

    /** Program entry point.

    • @param args program arguments (unused here)
      */
      public static void main(final String args) {
      try {
      // configure Orekit
      final File home = new File(System.getProperty(“user.home”));
      final File orekitData = new File(home, “orekit-data”);
      if (!orekitData.exists()) {
      System.err.format(Locale.US, “Failed to find %s folder%n”,
      orekitData.getAbsolutePath());
      System.err.format(Locale.US, “You need to download %s from %s, unzip it in %s and rename it ‘orekit-data’ for this tutorial to work%n”,
      “orekit-data-master.zip”, “https://gitlab.orekit.org/orekit/orekit-data/-/archive/master/orekit-data-master.zip”,
      home.getAbsolutePath());
      System.exit(1);
      }
      final DataProvidersManager manager = DataContext.getDefault().getDataProvidersManager();
      manager.addProvider(new DirectoryCrawler(orekitData));

       //  Initial state definition : date, orbit
 final AbsoluteDate initialDate = new AbsoluteDate(2023, 10, 26, 8, 00, 00.000, TimeScalesFactory.getUTC());
 final double mu =  3.986004415e+14; // gravitation coefficient
 final Frame inertialFrame = FramesFactory.getEME2000(); // inertial frame for orbit definition
 final Vector3D position  = new Vector3D(-7518565.475, 4279265.668, 4816446.308);
 final Vector3D velocity  = new Vector3D(-5463.904, -5533.9672, 50.719);
 final PVCoordinates pvCoordinates = new PVCoordinates(position, velocity);
 final Orbit initialOrbit = new KeplerianOrbit(pvCoordinates, inertialFrame, initialDate, mu);

 // Propagator : consider a simple Keplerian motion 
 final Propagator kepler = new KeplerianPropagator(initialOrbit);

 // Earth and frame
 final Frame earthFrame = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
 final BodyShape earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
                                              Constants.WGS84_EARTH_FLATTENING,
                                              earthFrame);

 //创建一个地面区域
 GeodeticPoint point1 = new GeodeticPoint(Math.toRadians(0), Math.toRadians(30), 0.0);
 GeodeticPoint point2 = new GeodeticPoint(Math.toRadians(30), Math.toRadians(30), 0.0);
 GeodeticPoint point3 = new GeodeticPoint(Math.toRadians(30), Math.toRadians(0), 0.0);
 GeodeticPoint point4 = new GeodeticPoint(Math.toRadians(0), Math.toRadians(0), 0.0);

 // 创建表示地面区域的SphericalPolygonsSet
 SphericalPolygonsSet region = EllipsoidTessellator.buildSimpleZone(1.0e-10, point1, point2, point3, point4);

 // Event definition
 final double maxcheck  = 60.0;
 final double threshold =  0.001;
 final double elevation = FastMath.toRadians(5.0);
  
 final GeographicZoneDetector zoneDetector = new GeographicZoneDetector(maxcheck, threshold, earth, region, elevation)
 		.withHandler((s, detector, increasing) -> {
             if (increasing) {
                 System.out.println("Visibility begins at " + s.getDate());
             } else {
                 System.out.println("Visibility ends at " + s.getDate());
             }
             return Action.CONTINUE; // Continue processing the event
         });

 // Add event to be detected
 kepler.addEventDetector(zoneDetector);

 // Propagate from the initial date to the first raising or for the fixed duration
 final SpacecraftState finalState = kepler.propagate(initialDate.shiftedBy(360000.));

 System.out.println(" Final state : " + finalState.getDate().durationFrom(initialDate));

      } catch (OrekitException oe) {
      System.err.println(oe.getLocalizedMessage());
      }
      }

}

18

This is what this detector computes.

19

The last parameter is not an elevation above ground, it is a margin counted on the Earth surface (i.e. similar to a longitude if East-West on equator or to a latitude if North-South).
What the detector computes is the signed distance from the satellite sub-point to the closest boundary point of the zone, and then it subtract the margin. With your 5° margin, you have roughly a zone that is 40° wide with rounded corners (30° + 5° on one side and 5° on the other side).

If you want the exact entry/exit of the zone, you should use 0° as the margin.

The margin is intended to be used when people need to be really sure they have fully entered the zone, in which case they will put a small negative margin like -0.01° to make sure they are really inside, or on the other hand when they want to lake sure they are away of a forbidden zone, then they would put a small positive value like +0.01°. In most cases, though, they want the zone itself so they put a 0° margin.

20

OK，thank you for your reply, I will give it a try.
Best wishes!

22

Why is there a roughly 40 ° wide rounded area for a 5 ° edge (30 °+5 ° on one side and 5 ° on the other side)?