Discrepencies during TargetPointingLaw

Orekit Python Wrapper usage
1

Hi everyone,

I’m trying to applicate switching conditions between sun pointing law and target Pointing law (TargetPointing)
The trigger is an elevation detector (Basically to detect a ground station visibility).
I’m then getting the OEM and AEM file (created with the handler, using CIC format) in order to put it back in VTS (3D visualization).

It worked well and I’m able to see the transition between sun pointing and Target pointing.
However, When the elevation is quite high, I can see a Discrepency between the targeted point and the true target point calculated by Orekit. So it miss the target (reasonably) as you can see in this picture.

pointing_gs
pointing_gs279×248 61.2 KB

I used of course the same frame in Orekit and VTS (EME 2000).
Here is my code :


ra = 500 * 1000.0  # Apogee in meters
rp = 500 * 1000.0  # Perigee in meters
# if SSO_Choice != 1, please indicate the orbit inclination, in degrees
i = 80  # Inclination in deg
omega = float(0)  # Perigee argument in deg
raan = float(0)  # Right ascension of ascending node in deg
lv = float(0)  # True anomaly in deg

a = (rp + ra + 2 * Constants.WGS84_EARTH_EQUATORIAL_RADIUS) / 2.0  # SMA
e = 1.0 - (rp + Constants.WGS84_EARTH_EQUATORIAL_RADIUS) / a  # Ecc

if SSO_Choice == 1:
    i = SSO_inclination_calculations(a,e)

# Satellite definition
satellite_mass = float(1.1)
area = float(0.015)

gs = GroundPoint("GS", latitude=45, longitude=20, altitude=10.0)

minimum_detectable_elev = 5.0

# ##################### END OF USER INPUTS #####################################
# ##############################################################################

propagator = PropagatorFactory.from_keplerian(
    a, e, i, omega, raan, lv, satellite_mass, area, initial_date
)

Elevation_Detector = Elevation_detector_wrapper(
    gs.topocentric, minimum_detectable_elev
)
Elevation_Detector.add_event_detectors_to_propagator(propagator)


sun = CelestialBodyFactory.getSun()
Sun_Pointing_law=CelestialBodyPointed(propagator.getFrame(), sun, Vector3D.MINUS_I, Vector3D.MINUS_K, Vector3D.MINUS_J)
GS_pointing_law = TargetPointing(propagator.getFrame(),gs.geodetic_pt,earth)

elevation_detec = (
            ElevationDetector(gs.topocentric)
            .withConstantElevation(FastMath.toRadians(minimum_detectable_elev))
            .withMaxCheck(10.0)
            .withThreshold(1.0e-6)
            .withHandler(ContinueOnEvent())
        )

time_step = 1.0
Handler = Handler_For_OEM_AEM(initial_date)     # The handler is used to get all parameters 
Propagator.cast_(propagator).setStepHandler(time_step,Handler)


attitudesSequence = AttitudesSequence()
attitudesSequence.addSwitchingCondition(Sun_Pointing_law,GS_pointing_law,
                                        elevation_detec, True, False, 0.10, 
                                        AngularDerivativesFilter.USE_RRA, None)
attitudesSequence.addSwitchingCondition(GS_pointing_law, Sun_Pointing_law, elevation_detec, 
                                        False, True, 0.10, AngularDerivativesFilter.USE_RRA,  None)


attitudesSequence.resetActiveProvider(Sun_Pointing_law)
propagator.setAttitudeProvider(attitudesSequence)

attitudesSequence.registerSwitchEvents(propagator)
# Propagate
final_state = propagator.propagate(
    initial_date, initial_date.shiftedBy(Propagation_duration))


OEM_creation_file(Handler)
AEM_creation_file(Handler)

And the Handler is defined by :


def Handler_For_OEM_AEM(first_date):
    class Prop_Handler(PythonOrekitFixedStepHandler):

        def init(self, so, t, step):
            super().__init__()
            self.data = {
                'dates': [],
                'dt': [],
                'Pos_X': [],
                'Pos_Y': [],
                'Pos_Z': [],
                'q0': [],
                'q1': [],
                'q2': [],
                'q3': [],
            }
            self.print_counter = 0

        def handleStep(self, state):
            date = state.getDate()
            x= state.getOrbit().getPosition(FramesFactory.getEME2000()).getX()
            y=state.getOrbit().getPosition(FramesFactory.getEME2000()).getY()
            z=state.getOrbit().getPosition(FramesFactory.getEME2000()).getZ()
            q0= state.getAttitude().getRotation().getQ0()
            q1=state.getAttitude().getRotation().getQ1()
            q2=state.getAttitude().getRotation().getQ2()
            q3=state.getAttitude().getRotation().getQ3()
            dt=date.durationFrom(first_date)

            if self.print_counter % 500 == 0:  # Print date every 500 steps (Saving calculations)
                print(date)
            mjd = date.getMJD()
            mjd_day = int(mjd)  # Partie entière
            mjd_sec = (mjd - mjd_day) * 86400  # Partie fractionnaire en secondes
            # Collect data and put them in a dataframe
            #ëdate_str = date.toString().split('.')[0]
            self.data['dates'].append(mjd_day)
            self.data['dt'].append(round(mjd_sec))
            self.data['Pos_X'].append(x/1000.0)
            self.data['Pos_Y'].append(y/1000.0)
            self.data['Pos_Z'].append(z/1000.0)
            self.data['q0'].append(q0)
            self.data['q1'].append(q1)
            self.data['q2'].append(q2)
            self.data['q3'].append(q3)
            self.print_counter += 1

        def finish(self, finalState):
            pass

Do you know what could cause a discrepency (for high elevation) ? Do you think it’s just related to VTS and not the Orekit calculation ?

Thanks and have a nice day !

Cheers,
Emilien

2

Good morning Emilien,

I noticed the same issue a while ago with files that were not generated with Orekit. It is indeed an issue with VTS, I checked with their support team.

The problem was that VTS used a spherical Earth and the oblateness was only taken into account in SurfaceView, not Celestia. A partial solution was to:

  • start the VTS simulation as normal
  • close celestia
  • edit the auto-generated Celestia file in Apps/Celestia/bin/extras_1/VTS/PROJECT.ssc. I can’t remember exactly but somewhere in this file there should be a setting for the oblateness of the Earth
  • in the Broker, Apps tab, unfold Celestia and uncheck include launcher start
  • restart Celestia from the Broker

However:

  • the lat/lon station coordinates correspond to a spherical Earth, so you will have to adjust them to use elliptical coordinates instead.
  • I think the textures are not perfectly aligned
  • I could still see a difference of a few hunders of meters in longitude, so there is probably also a discrepancy in the Earth rotation parameters.

This was probably with an old version of VTS, now looking at their changelog I see an entry for version 3.7.1 that suggests this was fixed, so maybe you are still using an older version ?

Clément

3

Hi clement,

Thank you for your answer. I’m using the latest version of VTS (Version 3.8.1). In the VTS documentation, the UseWGS84 reference system is described as follow : (you are right about the SurfaceView).

UseWGS84 reference system : in this mode, the flattening (ellipticity, or oblateness) of the body is taken into account. For the moment, only SurfaceView viewing Earth handles this configuration.

I tried your solution regarding the ssc file and found the following lines :

Modify "Earth" "Sol"
{
   Class "planet"
   OrbitColor [1 0 0.0601663]
   TrajectoryPlotWidth 1
   Radius 6378.137
   Oblateness 0.00335281066474748
}
  • the lat/lon station coordinates correspond to a spherical Earth, so you will have to adjust them to use elliptical coordinates instead.

Do you know how can I manage to do this ? I could then put these parameters in VTS to see if it changes.

4

Hi Émilien,

I think I used the CelestLab library to perform the conversion to elliptical coordinates, but I’m not sure.

In any case I would recommend that you reach out to the VTS support team because maybe some things have changed in recent versions.

Clément

5

Thank you for your help, @cmasson ! I’m going to clean up my code and post it back here to help other potential Orekit users.