THE PRECISE AUTONOMOUS ORBIT KEEPING EXPERIMENT ON THE PRISMA MISSION Sergio De Florio∗* and Simone D’Amico* This paper analyses the problem of the autonomous control of the Longitude of the Ascending Node (LAN) for a satellite in Low Earth Orbit (LEO) by means of along-track and anti-alongtrack velocity increments which adjust the semi-major axis. The problematic concerned with the possibility of generating the reference orbit (RO) on-board and with the estimation of the atmospheric drag are considered. The Autonomous Orbit Keeping (AOK) experiment of the PRISMA formation flying mission will be the test platform of the theoretical achievement here exposed. PRISMA comprises a fully maneuverable micro-satellite (MAIN) as well as a smaller sub-satellite (TARGET) which are launched together in a clamped configuration and separated in orbit after completion of all checkout operations. The AOK on-board software shall demonstrate MAIN satellite autonomous orbit control using a guidance law for the orbits LAN and shall implement a deterministic control algorithm using along-track and anti-along-track velocity increments. Using GPS-based absolute navigation data, AOK shall command thruster activations in the orbital frame to autonomously control the orbit within a predefined window. The AOK experiment represents a substantial contribution to the state of the art and paves the way to the accurate and autonomous orbit control of LEO satellites on a routine basis. The main requirement of the experiment is to demonstrate an orbit control accuracy of the osculating ascending node of 10 m (1σ). The paper shows results from real-world software simulations where the accuracy of the reference orbit is limited and GPS sensors and hydrazine actuators are accurately modeled. The fundamental approach on which the software design, validation and testing is based, is also explained.

INTRODUCTION Autonomous navigation and orbit control is acquiring increasing interest as it can enhance mission performance and provide significant operations cost reduction. By controlling the spacecraft orbit to match a chosen reference, the fulfillment of strict requirements on different orbit parameters can be achieved in real time and with a significant reduction of ground operations. The use of autonomous navigation and control can afford cost benefits by relieving of ground station operational burden, reduction of planning and scheduling costs, use of less propellant than traditional orbit maintenance, need of smaller and lighter weight thrusters, easier scheduling ground station operations and data collection (Ref. 1). The benefits of this type of control can be fully exploited in remote sensing space missions where the satellites are generally placed in sun-synchronous, phased and frozen orbits: examples of missions of this type are TerraSAR-X (Ref. 2, 3 and 4) and ERS-1, ERS-2 (Ref. 5 and 6). This paper analyses the problem of the autonomous control of the longitude of ascending node for a satellite in LEO by means of adjustments of the semi-major axis. Establishing a reference orbit and the way of generating it is a key element of the orbit control system. The on-board generation by a numerical routine can be severely constrained by the limitation of the on-board computer. On the other hand generating the reference points with an analytical model can be precluded by the accuracy requirements of the grid that defines the nominal longitudes of the ascending node during the satellite life. Uploading a reference orbit from the ground to the satellite keeps from facing all the problems concerned with on-board generation but reduces the autonomy and flexibility of the control system. As this type of orbit control is based on along-track and anti-along-track velocity increments which adjust the semi-major axis, a correct estimation of the atmospheric drag, the main non gravitational perturbation in LEO, is also essential to the end of a fine control. The autonomous orbit keeping experiment is the secondary objective of the German Aerospace Center (DLR) contributions to the PRISMA mission (Ref. 7 and 8). PRISMA is a micro-satellite mission created ∗ Space

Flight Technology Department, German Aerospace Center (DLR), D-82230 Wessling, Germany.

by the Swedish National Space Board (SNSB) and Swedish Space Corporation (SSC), which serves as a test platform for autonomous formation flying and rendezvous of spacecraft. PRISMA comprises a fully maneuverable micro-satellite (MAIN) as well as a smaller sub-satellite (TARGET) which are launched together in a clamped configuration and separated in orbit after completion of all checkout operations (Figure 1). The design orbit of the satellites is sun-synchronous, near-polar at approximately 700 km altitude with local time of the ascending node at 18.00. The mission schedule foresees a launch in 2009 of the two spacecraft with a targeted lifetime of at least eight months. The PRISMA mission primary objective is to demonstrate in-flight technology experiments related to autonomous formation flying, homing and rendezvous scenarios, precision close range 3D proximity operations, soft and smooth final approach and recede maneuvers, as well as to test instruments and unit developments related to formation flying. The backbone navigation sensor is based on GPS receivers on both satellites (Ref. 9). At the end of the mission, after the formation release, the AOK on-board software shall demonstrate MAIN satellite autonomous orbit control using a guidance law for the orbits LAN that implements an analytical feedback control algorithm using along-track and anti-along-track velocity increments. Using GPS-based absolute navigation data, AOK shall command thruster activations in the orbital frame to autonomously control the orbit within a predefined window. The main requirement of the experiment is to demonstrate a control accuracy of the osculating ascending node of 10 m (1σ). Main difference with respect to similar experiments (Ref. 10 and Ref. 11) are the extremely strict required control accuracy and the full autonomy enhanced by the possibility of on-board reference orbit propagation. The design of AOK, as that of the entire on-board PRISMA flight software, follows a model based approach. It makes use of Matlab/Simulink to describe the data flow, functionality and scheduling properties of individual software components. Following the conceptual design and the validation in a non real-time environment, a real-time executable of the flight software for the XPC platform (Ref. 12) or the LEON3 target processor (Ref. 13) is generated using the automatic code generator of Real Time Workshop. The full consistency of the flight software generation chain is verified by comparing relevant outputs obtained from application runs on a host standard laptop PC in a Matlab/Simulink environment and on the target LEON3 board under the Real-Time Executive for Multiprocessor Systems (RTEMS, Ref. 14).

Figure 1 PRISMA formation

THE PROBLEM OF LAN CONTROL In this section the general problem of LAN control is reviewed. The free motion of the real orbit with respect to the reference orbit is analyzed in order to identify which orbit perturbations are relevant in the scope of this specific control problem. The relation between the evolution in time of the LAN and the semimajor axis is then considered as this is the basis of the control algorithm that is to be implemented. Finally the problems of on-board atmospheric drag estimation and reference orbit generation are considered as these are key factors in the realization of a fine autonomous orbit control. Deviation of the Real Orbit from the Reference Table 1 shows the MAIN satellite nominal mean Keplerian elements. The mean orbital elements of a satellite in LEO deviate from their nominal values under the action of perturbing forces. Orbital elements values are dictated by specific mission requirements. In the case of MAIN, the orbit is required to be sun synchronous (in general orbits of remote sensing satellites, for which the AOK experiment is of special interest, are imposed to be sun-synchronous, phased and frozen at the same time). Table 1 MAIN satellite nominal orbital elements Keplerian Element

Value

Semi-major axis (a) Eccentricity (e) Inclination (i) Right Ascension of Ascending Node (Ω) Argument of perigee (ω) Mean anomaly (M)

7087297 (m) 0.0015 (-) 98.1877 (deg) 189.8914 (deg) 1.0938 (deg) -1.09369 (deg)

Sun synchronism is imposed through the relationship between a and i obtained by setting the secular angular rotation of the line of nodes, due to the non-sphericity of the Earth, equal to the known angular rotation of the meridian plane containing the mean Sun. The sun synchronism condition is thus represented by Eq. (1) where the secular variation of Ω can be represented with a good approximation by Eq. (2) that considers only the first three zonal terms J2 , J3 and J4 (Ref. 15): ˙ = ωSa Ω

˙ = − 3 nJ2 Ω 2

e2 4

µ

RE p

¶2

3 cos i − nJ22 2

µ

RE p

¶4

½ cos i

· ¸ 1 1 9 3 5 9 + (1 − e2 ) 2 − sin2 i + (1 − e2 ) 2 + 4 2 2 4

µ ¶ ¾ µ ¶3 5 e2 3 RE 1 + sin2 i + (7 − 15 sin2 i) cos 2ω − nJ3 (15 sin2 i − 4)e cot i sin ω + 4 8 8 p

15 nJ4 16

µ

RE p

¶4

(1)

(2)

¶ ¸ · µ 3 cos i (4 − 7 sin2 i) 1 + e2 − (3 − 7 sin2 i)e2 cos 2ω 2

where: ωSa = 1.99099299 × 10−7 rad s−1 is the mean apparent revolution speed of the Sun around the Earth

RE = 6378140 m is the equatorial radius of the terrestrial spheroid p = a(1 − e2 ) is the parameter of the orbit p n = µ/a3 is the mean motion (with µ = 3.9860064 × 1014 m3 s−2 ).

Figure 2 Deviation of LAN from reference orbit value

A reference orbit is an orbit representing the mean nominal motion of the satellite over a long time interval. Thus the reference orbit model, in the case here considered, should at least consider the non-spherical terms of the Earth gravity potential which cause the sun synchronous secular motion of the LAN. Other secular, long-periodic and short periodic terms can be included in the reference orbit depending on the flight control accuracy requirements. Deviations of real from reference orbital elements, leading to a violation of flight control requirements, have to be corrected by orbit maneuvers. Figures 2 and 3 show the computed contributions of the different perturbations to the deviation of the LAN (controlled parameter) and the semi-major axis (control parameter) from the reference orbit corresponding values. The analysis was performed by numerical orbit propagation over 35 days (the scheduled AOK experiment duration) and using orbit propagation model and satellite physical parameters given in Table 2. At the altitude of MAIN satellite (about 700 km) the major non-conservative perturbation force is the atmospheric drag whose main effect is the secular decrease of the semi-major axis and, as a consequence, the fast deviation of the LAN from its nominal value. The third body perturbation has a not negligible influence on the LAN deviation but a periodic and not appreciable effect on the semi-major axis. The solar radiation pressure (SRP) has a negligible influence on both the LAN and the semi-major axis. Relationship between LAN and Semi-major Axis Evolution A preliminary study of the relationship between the LAN time evolution i.e. the phase difference at the equator and the semi-major axis is fundamental for the orbit control strategy definition. The phase difference ∆LAN is defined as the difference between the real orbit track and the reference orbit track measured along the equator. Considering only the geopotential term J2 and linearizing around the nominal orbit values, the general analytical espression for the evolution of ∆LAN with the orbital elements is (Ref. 17):

Figure 3 Deviation of semi-major axis from reference orbit value

d2 3 (∆LAN ) = − (ωSi − Ω˙ R ) dt2 2

µ

RE a

¶"

# µ ¶2 ˙R 7 Ω 7 RE da 1+ + J (4 cos2 i − 1) ˙ R) 2 2 3 (ωSi − Ω a dt (3)

" −RE

µ ˙ R tan i + 6J2 Ω

RE a

¶2

#

˙ R ) sin 2i di (ωSi − Ω dt

where: ΩR is the right ascension of the ascending node of the reference orbit ωSi is the sidereal revolution speed of the Sun around the Earth. The angular speed ωSi is tied to the mean apparent revolution speed of the Sun around the Earth ωSa and to the angular rotation speed of the Earth ωE by the equation: ωSi = ωSa + ωE with ωSa = 2π/TSa ωE = 2π/TE = 7292115 × 10−11 rad s−1 . If the orbit is sun synchronous (Sy ) Ω˙ R = ωSa and considering Eq. (4), Eq. (3) becomes:

(4)

Table 2 Propagation parameters MAIN Physical Properties

Value

Mass Drag area Drag coefficient SRP effective area SRP coefficient

154.4 (kg) 1.3 (m2 ) 2.5 (-) 2.5 (m2 ) 1.3 (-)

Real Orbit Propagation

Value

Earth gravity feld Atmospheric density Sun and Moon ephemerides Solar Radiation Pressure

GRACE GGM01S 90x90 Harris-Priester Simplified Analytical Formulas (Ref. 16) Cannon-Ball Model

Reference Orbit Propagation

Value

Earth gravity field

GRACE GGM01S 90x90

d2 3π (∆LAN ) ≈ − dt2 TE

µ

RE a

¶

"

da di (1 + ²) − dt dt

µ

da di

¶

# (5)

Sy

with µ

²≈

da di 7 3

¶ Sy

µ

2 =− a 3

TE TSa

¶

µ

7 + J2 2 µ

η = 12J2

TSa TE

TE TSa µ

¶µ

RE a

¶µ

1+η 1+²

¶ tan i

¶2

RE a

(4 cos2 i − 1) ¶2 cos2 i

and η ¿ 1 and ² ¿1 for near polar orbits. Basic Orbit Keeping Strategy Referring to Eq. (3), if da/dt and di/dt are taken as constants, ∆LAN is found to have a near parabolic variation as its second derivative is constant. This means that if it is assumed that the semi-major axis has a linear decay under the influence of the atmospheric drag, the real orbit LAN will move parabolically eastwards of the reference orbit LAN. Based on this consideration a simple maneuver cycle can be defined for the ideal case: when the deviation between the real and nominal LAN exceeds a pre-defined upper bound ∆LM AX a correcting impulse in the along-track direction of the satellite’s orbit is applied as much as twice as necessary for recovering the initial semi-major axis value. As a result, just after the impulse, the satellite’s LAN begins to move westward, reaches the lower bound of the allowable band for the LAN and drifts back to the upper limit, where the next correction maneuver is made. The theoretical maneuver cycle and the semi-major increments to be applied can thus be obtained from Eq. (5) assuming also that, as the reference orbit is sun-synchronous, di/dt can be taken as negligible compared to da/dt.

s Tc =

µ Kc

a RE

s

¶ ¯ ¯−1 ¯ da ¯ ¯ ¯ ¯ dt ¯

∆ac =

µ Kc

a RE

¶¯ ¯ ¯ da ¯ ¯ ¯ ¯ dt ¯

(6)

with Kc = (16/3π)TE ∆LM AX

Figure 4 Semi-major axis deviation (top), LAN deviation (middle), magnitude of along-track maneuvers (bottom) computed by AOK during the guidance process. Test scenario: constant drag, no other perturbation

The velocity increment to be applied is: µ ∆Vc =

∆ac 2a

¶ V

(7)

where V is the mean velocity. Figure 4 shows an example of guidance process in the ideal case in which the atmospheric density has been fixed at the value 12 · 10−14 kg/m3 (the average of the Harris-Priester model at 700 km altitude) and the atmospheric drag is the only orbit perturbation. At each maneuver cycle of about 1.3 days the evolution of the LAN deviation is evenly parabolic while that of the semi-major is linear. The simulation is performed over 35 days. Figure 5 shows the guidance process in a test scenario in which still the atmospheric density is held constant but Sun-Moon, radiation pressure and tides perturbations are considered. In this case some unevenness is introduced mainly by the effect of the Sun-Moon perturbation action. Atmospheric Drag Estimation By considering Equation (6) it is straightforward to conclude that a correct and fine estimation of the correction maneuver is strictly connected to a correct estimation of da/dt which, as already outlined (Figure 3),

Figure 5 Semi-major axis deviation (top), LAN deviation (middle), magnitude of along-track maneuvers (bottom) computed by AOK during the guidance process. Test scenario: constant drag, Sun-Moon, solar radiation pressure, tides and relativity gravity field perturbations

is mainly determined by the atmospheric drag. As at the altitude of MAIN the atmospheric density has a remarkable variability with time, mainly due to the variability of the solar flux, an atmospheric model cannot give the required accuracy. On the other hand the accelerometers on-board the MAIN satellite are not suitable for the estimation of the atmospheric drag as the typical bias is 4 · 10−2 m/s2 , much larger than typical atmospheric drag acceleration values at 700 km (< 10−5 m/s2 ). The problem has been solved starting by the simple consideration that all the information required for the estimation of the da/dt is available on-board the spacecraft: the knowledge of the actual orbit estimated by the on-board navigation filter with an accuracy of 2 m 1σ (Ref. 9) and the reference orbit (uploaded or computed on-board). The value of da/dt is then found with the following steps: 1. Filter the actual and reference orbit semi-major axis data by taking their values at the ascending node passes. 2. Compute the differences ∆aAN between the actual and reference orbit semi-major axis values. 3. Smoothing of data if required. 4. Find da/dt as the angular coefficient of the line representing the linear fit over the time of a sample of ∆aAN data. Figure 6 shows this procedure applied in the case in which the real orbit is generated by numerical orbit propagation including the aspheric Earth gravity field (GGM01S model) through an expansion in spherical harmonics up to degree and order 90, tides and relativity gravity field perturbations, the Sun and Moon third body forces, the atmospheric drag (atmospheric density model Harris-Priester) and solar radiation pressure. The reference trajectory has been generated using the same gravity field model used to generate the real orbit. Navigation errors are not included. Figure 7 shows the corresponding guidance process in which a sample buffer of 10 ∆aAN ’s is renewed at each ascending node passage for the estimation of da/dt. The buffer is

Figure 6 ∆aAN data fitting process

reset after each maneuver. It can be noticed that the first maneuver is over-estimated because at this point only 5 samples were available and, as no fitting was possible, the da/dt has been estimated using a backup atmospheric model.

Reference Orbit Generation The reference orbit can be generated either on-ground and uploaded on-board the satellite during a ground station pass or directly on-board. The great advantage of the on-ground generation is that a numerical orbit propagation including all the available harmonics of the gravity field model can be used thus achieving the best control performance. The disadvantage is that the reference orbit has to be uploaded periodically (3 days in the case of MAIN) imposing a short time limit on the autonomy of the orbit control. On the other end while an on-board reference orbit generation allows much longer autonomous control periods, the performance is degraded by the fact that only limited gravity field models can be exploited due to the limitation imposed by available on-board CPU resources. In the case that a poor on-board generated orbit is used, the estimation of da/dt turns to be difficult. The short period oscillation of the ∆aAN cannot be catched due to the lack of short period perturbations in the gravity field model. Figure 8 shows the comparison of the ∆aAN curves obtained with a reference orbit generated with a 90x90 gravity field model and one generated with a 5x5 model. In the case of the 5x5 orbit model a short period (for example 10 points) linear fitting is totally unreliable for estimating the atmospheric drag. A long period (50 points) linear fitting would allow an accurate estimation of da/dt, but is not useful at the end of the orbit control because it would require a minimum maneuver cycle period of 50 orbits (i.e. more than 3 days). Nevertheless an efficient estimation of the da/dt is still possible using samples of 10 ∆aAN values but this time applying a proper smoothing process. The smooth filter takes as input the points marked by circles in Figure 8 and flattens them on a line λ keeping unaltered the ratios of their distances from it. Line λ is defined by point (x1 , y) and angular coefficient m where x1 is the x-coordinate of the first point, y is the mean of the y-coordinates of all the points of the sample and m is the estimated da/dt. A guess value of da/dt is used for the very first smoothing process. The output of the smooth filter are the points marked by triangles.

Figure 7 Semi-major axis deviation (top), LAN deviation (middle), magnitude of along-track maneuvers (bottom) computed by AOK during the guidance process (external uploaded reference orbit). Test scenario: drag (Harris-Priester model), Sun-Moon, solar radiation pressure, tides and relativity gravity field perturbations

AUTONOMOUS ORBIT KEEPING The AOK experiment is a minimal target control problem as it has to demonstrate the autonomous fine control of the LAN with respect to a predefined reference orbit. The control accuracy with respect to a predefined phasing grid fixed to a given reference frame on Earth (for example a given longitude) can then be assessed by analyzing the relative movement of the reference orbit and the phasing grid. AOK operations will commence after the end of the PRISMA nominal operations. The test campaign foresees 1 week of initial configuration, calibration and open-loop tests and 4 weeks of close-loop operations. In this section the AOK controller logic and the software development and testing environment are described. AOK Features The AOK control cycle is computed from Equations (6) with da/dt estimated through the process described in the preceding section. The main control logic is the following: 1. An along-track velocity increment ∆V ( Eq. (7) ) is commanded if the time elapsed from last maneuver is larger than a pre-defined time interval ∆tV + , the LAN deviation exceed a pre-defined upper limit (∆LAN > ∆LM AX ) and ∆LAN is growing with time (d∆LAN /dt > 0). Constants ∆LM AX and ∆tV + are uploaded to the spacecraft via Telecommand (TC). 2. After the first along-track maneuver in order to calibrate ∆ac on the actual atmospheric drag value, it is recomputed using Eq. (8): ¯ ¯ ¯ da ¯ Tac ∆ac ∆acal = ¯¯ ¯¯ + dt ac 2 2

(8)

Figure 8 Estimation of da/dt for a reference orbit generated with a 5x5 gravity field model

where: (da/dt)ac is estimated by the long period fitting (with the number of points imposed by TC) if possible or by the ordinary short period fitting, and ∆acal is the calibrated value of ∆ac . Tac is the actual cycle period i.e. the time elapsed from the actuation of last maneuver. 3. An anti-along-track velocity increment is commanded if this possibility is enabled by TC, the LAN deviation exceeds a pre-defined lower limit (∆LAN < −n·∆LM AX where n is a real constant), d∆LAN /dt < 0 and the time elapsed from last maneuver is larger than a pre-defined time interval ∆tV − . Constants n and ∆tV − are uploaded to the spacecraft via TC. In this case: ∆acal

¯ ¯ ¯ da ¯ = ¯¯ ¯¯ dt

ac

Tc B

(9)

where: Tc is computed with Eq. (6) B is an empirical correction factor imposed by TC. Figure 7, provides an example of the AOK guidance process in case of a reference orbit is uploaded from ground. ∆LM AX is 10 m, the anti-along-track maneuvers are not enabled, ∆tV + is zero. As the reference orbit is the best available and neither navigation nor sensors and actuators errors are introduced in the simulation, ∆LAN has a smooth parabolic evolution. The mean of ∆LAN is -5 m and σ (standard deviation) is 12.3 m: the control accuracy of the osculating ascending node is then compliant with the main requirement (10 m 1σ). In contrast Figure 9 gives an example of AOK guidance with an on-board propagated orbit (gravitational field model GGM01S 5x5). ∆LM AX is 50 m and ∆tV + is 5 orbital periods (about 8 hours). The enlarged control window and waiting times are chosen to compensate the presence of oscillations in the evolution of ∆LAN which are due to the incompleteness of the reference orbit. The enlargement of this two parameters allows AOK to correctly estimate da/dt. A parabolic but noisy evolution of ∆LAN can still

Figure 9 Semi-major axis deviation (top), LAN deviation (middle), magnitude of along-track maneuvers (bottom) computed by AOK during the guidance process: on-board generated reference orbit, ∆LMAX = 50 m, ∆tV+ = 0

Figure 10 Semi-major axis deviation (top), LAN deviation (middle), magnitude of along-track maneuvers (bottom) computed by AOK during the guidance process: on-board generated reference orbit, ∆LMAX = +50 m ÷ -100 m, ∆tV+ = 8 h, ∆tV− = 0

be observed. This time the mean of ∆LAN is -50 m and σ is 53 m: the control accuracy is not compliant with the main requirement. Figure 10 shows the same test case as Figure 9 but with the enabling of anti-along-track maneuvers: the control window is this time bounded between -100 m and +50 m ( ∆LM AX = 50 m, n = 2). In this case the performance of AOK is better (the mean of of ∆LAN is -22 m and σ is 45 m) but at the cost of a larger number of maneuvers which impose a negative ∆V most of the times. Software Development and Testing Environment Figure 11 shows the simplified architecture of the PRISMA on-board software (OBS) architecture for the MAIN satellite and the simulation environment (see also Ref. 18). The PRISMA OBS has a layered structure with a Basic Software (BSW) level and an Application Software (ASW) level communicating with each other through dedicated message queues. While the BSW includes basic applications, device drivers and I/Outilities, the ASW encapsulates all top-level applications like spacecraft navigation, control, telecommand and telemetry. The GPS-based navigation system is split into three modules located in different OBS levels and running at different sample rates. The GPS interface (GIF) is part of the BSW, runs at 1 s sample time and is directly fed with GPS messages issued by the Phoenix GPS receivers emulator (PEM). GIF handles GPS raw data formats and ephemerides, and performs data sampling as well as coarse editing prior to the GPS-based orbit determination. The GPS-based Orbit Determination (GOD) and GPS-based Orbit Prediction (GOP) are embedded in the ASW layer as part of the ORB core (30 s sample time) and the GNC core (1 s sample time), respectively. GOD implements an extended Kalman filter to process GRAPHIC observables as well as single difference carrier phase measurements from MAIN satellite. Attitude data from both spacecraft are applied to correct for the GPS receivers antenna offset with respect to the spacecraft center of mass. Furthermore, a history of maneuver data is provided to GOD and taken into account in the orbit determination task. GOD performs a numerical orbit propagation which is invoked after the measurement update and provides orbit coefficients for interpolation to GOP for both spacecraft. The GOP module interpolates the orbit coefficients provided by GOD and finally supplies the various GNC core functions as well as the PRISMA payload with continuous position and velocity data. Due to the different data rates of the GPS-based navigation modules, orbit maneuver data have to be taken into account in both GOD and GOP (cf. Figure 11). In particular at each GNC step, the GOP task accounts for maneuvers which have not been considered by GOD in the last orbit determination/prediction process. The AOK and Autonomous Formation Control (AFC) modules are located in the ORB and GNC cores respectively. Finally the orbit propagator (PROP) generates the MAIN simulated orbit. For PRISMA, a model-based design method is used for the complete on-board application software. While traditional software development is performed on a programming level and focused on software engineering, the model-based approach emphasizes systems engineering methods.

Figure 11 Simplified scheme of the simulation software architecture

This raises the abstraction level for the development in order to allow an efficient handling of ever more complex systems. Benefits of a model-based over a traditional development approach are expected in the field of robustness, predictability, scalability, efficiency and faster development schedules. The software is implemented in Matlab/Simulink as well as C/C++. In particular, C++ low level routines, which are interfaced with Matlab/Simulink high level structures, represents the computational layer of the software system including for example numerical integrations, data filtering, etc. Matlab/Simulink provides instead the communication layer, including for example input/output interfaces, complex logics, time synchronization. The prototyping of the software as well as the related analysis and simulations are performed first on a host standard laptop PC in a Matlab/Simulink environment. All the functions can then be auto-coded using the Matlab tool Real Time Workshop and executed under the operating system Real Time Executive for Multiprocessor Systems (RTEMS) on the target LEON3 processor (the PRISMA satellites on board computer). The full consistency of the flight software generation chain can thus been verified by comparing relevant outputs generated by the software running on the host and target computer. One of the greatest advantages of this straightforward and cheap approach is that it allows the test of the flight software performance directly on the on-board spacecraft computer. Table 3 AOK execution times and CPU loads on LEON3 board AOK Mode

Execution Time (ms)

CPU Load (%)

63 596 49623

0.2 2 165

External uploaded RO On-board generated RO (multi-step) On-board generated RO (one-step)

The analysis of the software execution times has special interest because it allows to verify if the different software modules comply with the maximal allocated CPU loads. Analysis of heap and stack allocation and software profiling can also be made for the detection of eventual software execution bottlenecks. Table 3 collates the execution times and CPU loads of AOK running on the LEON3 board in three modes: with an external uploaded RO, a RO generated on-board by a propagation performed in one call or in several calls. Only the one-step propagation is not compliant with the CPU load constraints. External uploaded and multi-step on-board propagated RO modes are fully compliant with the CPU load constraints because AOK is located in the ORB core that is a low priority task. The possibility of using AOK in both external and internal propagated RO modes, enhances the robustness and the autonomy of the control as AOK can automatically switch to on-board generated RO in case of not validity of the external uploaded RO. Numerical Simulation Figure 12 shows how the guidance process of AOK in the test scenario of Figure 7 changes if the navigation errors are included and sensors and actuators are modeled by applying realistic errors to measurements and control actions respectively. Table 4 collects the navigation accuracy and the parameters used for the accuracy modeling of sensors and actuators. The propulsion system is characterized by a Minimum Impulse Value (MIV) of 0.7 mm/s and a Minimum Impulse Bit (MIB) of 0.07 mm/s. Consequently, the thrusters can only realize ∆V ’s which are larger than MIV and integer multiples of MIB. Furthermore, the performance error of the thruster system is quantified by Equation (10) ¯ ¯ ¯ ∆Vr − ∆Vp ¯ ¯ ¯ · 100 ²=¯ ¯ ∆Vp

(10)

where ∆Vp is the planned velocity increment computed by the AOK controller and ∆Vr is the actual velocity increment realized by the propulsion system. Finally the attitude control error, which causes deviations from the desired thrust direction, is treated as Gaussian noise with zero bias and 0.2◦ standard deviation in three axes. The following conclusions follow from an analysis of the simulation results shown in Figure 12:

Figure 12 Semi-major axis deviation (top), LAN deviation (middle), magnitude of along-track maneuvers (bottom) computed by AOK during the guidance process

1. The inclusion of navigation errors causes an oscillating evolution of ∆a and consequently a different evaluation of the LAN deviation with respect to the case in which no errors are included. It follows that the maneuver planning is different in the two cases. 2. The inclusion of navigation errors decreases the estimation accuracy of da/dt. It follows that, as the efficiency of the commanded maneuvers is reduced, a larger number of thrust pulses is needed in the guidance process. 3. The inclusion of sensors and actuators errors further reduces the efficiency of the commanded maneuvers. 4. In the case in which all the errors are included the mean of ∆LAN is -5 m and σ is 17 m. Table 4 Navigation, sensors and actuators accuracy modeling Absolute Navigation Accuracy

Value

Mean σ

0 (m) 2 (m)

Attitude Control Accuracy

Value

Mean σ

0 (◦ ) 0.2 (◦ )

Propulsion System Accuracy MIV MIB ²

Value −4

7 · 10 (m/s) 7 · 10−5 (m/s) 10

CONCLUSION The AOK experiment is the secondary objective of the DLR contributions to the PRISMA formation flying mission. It has to demonstrate autonomous orbit control of the longitude of the ascending node of a single spacecraft with an accuracy of 10 m (1σ). The software implements an analytical feedback control algorithm using along-track and anti-along-track velocity increments. In the realization of the AOK controller the problem of a correct estimation of the atmospheric drag has been solved. The software applies a filtering, linear fitting and smoothing of the difference between the reference and real orbits’ semi-major axis values at the ascending node passes. Advantages and disadvantages of using an on-ground or on-board generated reference orbit have been analyzed in terms of performance and operativeness. The results of the simulation performed are encouraging in showing that the fulfillment of the strict control accuracy requirement is feasible. Execution time tests on a representative PRISMA on-board computer indicate that at the present state the software is fully compliant with the CPU load requirements.

REFERENCES [1] J. R. Wertz, “Autonomous Navigation and Autonomous Orbit Control in Planetary Orbits as a Means of Reducing Operations Cost,” 5th International Symposium on Reducing the Cost of Spacecraft Ground Systems and Operations, Pasadena, CA USA, July 2003. [2] S. Kudryavtsev and E. Gill, “Autonomous Orbit Control of Remote Sensing Spacecraft using Spaceborne GPS Receivers,” 14 December 2000. [3] S. D’Amico, C. Arbinger, M. Kirschner, and S. Campagnola, “Generation of an Optimum Target Trajectory for the TerraSAR-X Repeat Observation Satellite,” 18th International Symposium on Space Flight Dynamics, Munich, Germany, 11-15 October 2004. [4] C. Arbinger, S. D’Amico, and M. Eineder, “Precise Ground-In-the-Loop Orbit Control for Low Earth Observation Satellites,” 18th International Symposium on Space Flight Dynamics, Munich, Germany, 11-15 October 2004. [5] M. Rosengren, “The Orbit Control of ERS-1 and ERS-2 for a very accurate Tandem Configuration,” RBCM, Vol. XXI Special Issue, 1999, pp. 72–78. [6] M. Rosengren, “Improved Technique for Passive Eccentricity Control,” AAS 89-155, AAS Symposium on Mission Design and Orbital Mechanics, 1989. [7] S. Persson, P. Bodin, E. Gill, J. Har, and J. J¨orgensen, “Prisma An Autonomous Formation Flying Mission,” Sardinia, Italy, ESA Small Satellite Systems and Services Symposium (4S), 25-29 September 2006. [8] E. Gill, S. D’Amico, and O. Montenbruck, “Autonomous Formation Flying for the PRISMA Mission,” Journal of Spacecraft and Rockets, Vol. 44, May 2007, pp. 671–681. [9] S. D’Amico, E. Gill, M. Garcia, and O. Montenbruck, “GPS-Based Real-Time Navigation for the PRISMA Formation Flying Mission,” Noordwijk, 3rd ESA Workshop on Satellite Navigation User Equipment Technologies - NAVITEC, 11-13 December 2006. [10] A. Lamy, M.-C. Charmeau, D. Laurichesse, M. Grondin, and R. Bertrand, “Experiment of autonomous Orbit Control on the Demeter Satellite,” 2004. [11] A. Lamy, M.-C. Charmeau, and D. Laurichesse, “Experiment of autonomous Orbit Control on the Demeter Satellite: In-flight Results and Perspectives,” AIAA Astrodynamics Specialist Conference, Keystone, Colorado, USA, 21-24 August 2006. [12] xPC Target 3 - Users Guide. The MathWorks, March 2008. [13] J. Gaisler, The LEON Processor User’s Manual V2.3.7. Gaisler Research, 2001. [14] RTEMS C Users Guide. Edition 4.6.5, for RTEMS 4.6.5. On-Line Applications Research Corporation (OAR) http://www.rtems.com, 30 August 2003. [15] J. Liu, “Satellite Motion about an Oblate Earth,” AIAA Journal, Vol. 12, November 1974, pp. 1511–1516. [16] O. Montenbruck and E. Gill, Satellite Orbits - Model, Methods and Applications. Heidelberg: Springer Verlag, 2000. [17] P. Micheau, Spaceflight Dynamics, Vol. I, ch. Orbit Control Techniques for LEO Satellites. Toulouse, France: Cepadues-Editions, 1995. [18] S. D’Amico, S. D. Florio, and T. Yamamoto, “Offline and Hardware-in-the-loop Validation of the GPS-based Real-Time Navigation System for the PRISMA Formation Flying Mission,” 3rd International Symposium on Formation Flying, Missions and Technology, ESA/ESTEC, Noordwijk, The Netherlands, 23-25 April 2008.

INTRODUCTION Autonomous navigation and orbit control is acquiring increasing interest as it can enhance mission performance and provide significant operations cost reduction. By controlling the spacecraft orbit to match a chosen reference, the fulfillment of strict requirements on different orbit parameters can be achieved in real time and with a significant reduction of ground operations. The use of autonomous navigation and control can afford cost benefits by relieving of ground station operational burden, reduction of planning and scheduling costs, use of less propellant than traditional orbit maintenance, need of smaller and lighter weight thrusters, easier scheduling ground station operations and data collection (Ref. 1). The benefits of this type of control can be fully exploited in remote sensing space missions where the satellites are generally placed in sun-synchronous, phased and frozen orbits: examples of missions of this type are TerraSAR-X (Ref. 2, 3 and 4) and ERS-1, ERS-2 (Ref. 5 and 6). This paper analyses the problem of the autonomous control of the longitude of ascending node for a satellite in LEO by means of adjustments of the semi-major axis. Establishing a reference orbit and the way of generating it is a key element of the orbit control system. The on-board generation by a numerical routine can be severely constrained by the limitation of the on-board computer. On the other hand generating the reference points with an analytical model can be precluded by the accuracy requirements of the grid that defines the nominal longitudes of the ascending node during the satellite life. Uploading a reference orbit from the ground to the satellite keeps from facing all the problems concerned with on-board generation but reduces the autonomy and flexibility of the control system. As this type of orbit control is based on along-track and anti-along-track velocity increments which adjust the semi-major axis, a correct estimation of the atmospheric drag, the main non gravitational perturbation in LEO, is also essential to the end of a fine control. The autonomous orbit keeping experiment is the secondary objective of the German Aerospace Center (DLR) contributions to the PRISMA mission (Ref. 7 and 8). PRISMA is a micro-satellite mission created ∗ Space

Flight Technology Department, German Aerospace Center (DLR), D-82230 Wessling, Germany.

by the Swedish National Space Board (SNSB) and Swedish Space Corporation (SSC), which serves as a test platform for autonomous formation flying and rendezvous of spacecraft. PRISMA comprises a fully maneuverable micro-satellite (MAIN) as well as a smaller sub-satellite (TARGET) which are launched together in a clamped configuration and separated in orbit after completion of all checkout operations (Figure 1). The design orbit of the satellites is sun-synchronous, near-polar at approximately 700 km altitude with local time of the ascending node at 18.00. The mission schedule foresees a launch in 2009 of the two spacecraft with a targeted lifetime of at least eight months. The PRISMA mission primary objective is to demonstrate in-flight technology experiments related to autonomous formation flying, homing and rendezvous scenarios, precision close range 3D proximity operations, soft and smooth final approach and recede maneuvers, as well as to test instruments and unit developments related to formation flying. The backbone navigation sensor is based on GPS receivers on both satellites (Ref. 9). At the end of the mission, after the formation release, the AOK on-board software shall demonstrate MAIN satellite autonomous orbit control using a guidance law for the orbits LAN that implements an analytical feedback control algorithm using along-track and anti-along-track velocity increments. Using GPS-based absolute navigation data, AOK shall command thruster activations in the orbital frame to autonomously control the orbit within a predefined window. The main requirement of the experiment is to demonstrate a control accuracy of the osculating ascending node of 10 m (1σ). Main difference with respect to similar experiments (Ref. 10 and Ref. 11) are the extremely strict required control accuracy and the full autonomy enhanced by the possibility of on-board reference orbit propagation. The design of AOK, as that of the entire on-board PRISMA flight software, follows a model based approach. It makes use of Matlab/Simulink to describe the data flow, functionality and scheduling properties of individual software components. Following the conceptual design and the validation in a non real-time environment, a real-time executable of the flight software for the XPC platform (Ref. 12) or the LEON3 target processor (Ref. 13) is generated using the automatic code generator of Real Time Workshop. The full consistency of the flight software generation chain is verified by comparing relevant outputs obtained from application runs on a host standard laptop PC in a Matlab/Simulink environment and on the target LEON3 board under the Real-Time Executive for Multiprocessor Systems (RTEMS, Ref. 14).

Figure 1 PRISMA formation

THE PROBLEM OF LAN CONTROL In this section the general problem of LAN control is reviewed. The free motion of the real orbit with respect to the reference orbit is analyzed in order to identify which orbit perturbations are relevant in the scope of this specific control problem. The relation between the evolution in time of the LAN and the semimajor axis is then considered as this is the basis of the control algorithm that is to be implemented. Finally the problems of on-board atmospheric drag estimation and reference orbit generation are considered as these are key factors in the realization of a fine autonomous orbit control. Deviation of the Real Orbit from the Reference Table 1 shows the MAIN satellite nominal mean Keplerian elements. The mean orbital elements of a satellite in LEO deviate from their nominal values under the action of perturbing forces. Orbital elements values are dictated by specific mission requirements. In the case of MAIN, the orbit is required to be sun synchronous (in general orbits of remote sensing satellites, for which the AOK experiment is of special interest, are imposed to be sun-synchronous, phased and frozen at the same time). Table 1 MAIN satellite nominal orbital elements Keplerian Element

Value

Semi-major axis (a) Eccentricity (e) Inclination (i) Right Ascension of Ascending Node (Ω) Argument of perigee (ω) Mean anomaly (M)

7087297 (m) 0.0015 (-) 98.1877 (deg) 189.8914 (deg) 1.0938 (deg) -1.09369 (deg)

Sun synchronism is imposed through the relationship between a and i obtained by setting the secular angular rotation of the line of nodes, due to the non-sphericity of the Earth, equal to the known angular rotation of the meridian plane containing the mean Sun. The sun synchronism condition is thus represented by Eq. (1) where the secular variation of Ω can be represented with a good approximation by Eq. (2) that considers only the first three zonal terms J2 , J3 and J4 (Ref. 15): ˙ = ωSa Ω

˙ = − 3 nJ2 Ω 2

e2 4

µ

RE p

¶2

3 cos i − nJ22 2

µ

RE p

¶4

½ cos i

· ¸ 1 1 9 3 5 9 + (1 − e2 ) 2 − sin2 i + (1 − e2 ) 2 + 4 2 2 4

µ ¶ ¾ µ ¶3 5 e2 3 RE 1 + sin2 i + (7 − 15 sin2 i) cos 2ω − nJ3 (15 sin2 i − 4)e cot i sin ω + 4 8 8 p

15 nJ4 16

µ

RE p

¶4

(1)

(2)

¶ ¸ · µ 3 cos i (4 − 7 sin2 i) 1 + e2 − (3 − 7 sin2 i)e2 cos 2ω 2

where: ωSa = 1.99099299 × 10−7 rad s−1 is the mean apparent revolution speed of the Sun around the Earth

RE = 6378140 m is the equatorial radius of the terrestrial spheroid p = a(1 − e2 ) is the parameter of the orbit p n = µ/a3 is the mean motion (with µ = 3.9860064 × 1014 m3 s−2 ).

Figure 2 Deviation of LAN from reference orbit value

A reference orbit is an orbit representing the mean nominal motion of the satellite over a long time interval. Thus the reference orbit model, in the case here considered, should at least consider the non-spherical terms of the Earth gravity potential which cause the sun synchronous secular motion of the LAN. Other secular, long-periodic and short periodic terms can be included in the reference orbit depending on the flight control accuracy requirements. Deviations of real from reference orbital elements, leading to a violation of flight control requirements, have to be corrected by orbit maneuvers. Figures 2 and 3 show the computed contributions of the different perturbations to the deviation of the LAN (controlled parameter) and the semi-major axis (control parameter) from the reference orbit corresponding values. The analysis was performed by numerical orbit propagation over 35 days (the scheduled AOK experiment duration) and using orbit propagation model and satellite physical parameters given in Table 2. At the altitude of MAIN satellite (about 700 km) the major non-conservative perturbation force is the atmospheric drag whose main effect is the secular decrease of the semi-major axis and, as a consequence, the fast deviation of the LAN from its nominal value. The third body perturbation has a not negligible influence on the LAN deviation but a periodic and not appreciable effect on the semi-major axis. The solar radiation pressure (SRP) has a negligible influence on both the LAN and the semi-major axis. Relationship between LAN and Semi-major Axis Evolution A preliminary study of the relationship between the LAN time evolution i.e. the phase difference at the equator and the semi-major axis is fundamental for the orbit control strategy definition. The phase difference ∆LAN is defined as the difference between the real orbit track and the reference orbit track measured along the equator. Considering only the geopotential term J2 and linearizing around the nominal orbit values, the general analytical espression for the evolution of ∆LAN with the orbital elements is (Ref. 17):

Figure 3 Deviation of semi-major axis from reference orbit value

d2 3 (∆LAN ) = − (ωSi − Ω˙ R ) dt2 2

µ

RE a

¶"

# µ ¶2 ˙R 7 Ω 7 RE da 1+ + J (4 cos2 i − 1) ˙ R) 2 2 3 (ωSi − Ω a dt (3)

" −RE

µ ˙ R tan i + 6J2 Ω

RE a

¶2

#

˙ R ) sin 2i di (ωSi − Ω dt

where: ΩR is the right ascension of the ascending node of the reference orbit ωSi is the sidereal revolution speed of the Sun around the Earth. The angular speed ωSi is tied to the mean apparent revolution speed of the Sun around the Earth ωSa and to the angular rotation speed of the Earth ωE by the equation: ωSi = ωSa + ωE with ωSa = 2π/TSa ωE = 2π/TE = 7292115 × 10−11 rad s−1 . If the orbit is sun synchronous (Sy ) Ω˙ R = ωSa and considering Eq. (4), Eq. (3) becomes:

(4)

Table 2 Propagation parameters MAIN Physical Properties

Value

Mass Drag area Drag coefficient SRP effective area SRP coefficient

154.4 (kg) 1.3 (m2 ) 2.5 (-) 2.5 (m2 ) 1.3 (-)

Real Orbit Propagation

Value

Earth gravity feld Atmospheric density Sun and Moon ephemerides Solar Radiation Pressure

GRACE GGM01S 90x90 Harris-Priester Simplified Analytical Formulas (Ref. 16) Cannon-Ball Model

Reference Orbit Propagation

Value

Earth gravity field

GRACE GGM01S 90x90

d2 3π (∆LAN ) ≈ − dt2 TE

µ

RE a

¶

"

da di (1 + ²) − dt dt

µ

da di

¶

# (5)

Sy

with µ

²≈

da di 7 3

¶ Sy

µ

2 =− a 3

TE TSa

¶

µ

7 + J2 2 µ

η = 12J2

TSa TE

TE TSa µ

¶µ

RE a

¶µ

1+η 1+²

¶ tan i

¶2

RE a

(4 cos2 i − 1) ¶2 cos2 i

and η ¿ 1 and ² ¿1 for near polar orbits. Basic Orbit Keeping Strategy Referring to Eq. (3), if da/dt and di/dt are taken as constants, ∆LAN is found to have a near parabolic variation as its second derivative is constant. This means that if it is assumed that the semi-major axis has a linear decay under the influence of the atmospheric drag, the real orbit LAN will move parabolically eastwards of the reference orbit LAN. Based on this consideration a simple maneuver cycle can be defined for the ideal case: when the deviation between the real and nominal LAN exceeds a pre-defined upper bound ∆LM AX a correcting impulse in the along-track direction of the satellite’s orbit is applied as much as twice as necessary for recovering the initial semi-major axis value. As a result, just after the impulse, the satellite’s LAN begins to move westward, reaches the lower bound of the allowable band for the LAN and drifts back to the upper limit, where the next correction maneuver is made. The theoretical maneuver cycle and the semi-major increments to be applied can thus be obtained from Eq. (5) assuming also that, as the reference orbit is sun-synchronous, di/dt can be taken as negligible compared to da/dt.

s Tc =

µ Kc

a RE

s

¶ ¯ ¯−1 ¯ da ¯ ¯ ¯ ¯ dt ¯

∆ac =

µ Kc

a RE

¶¯ ¯ ¯ da ¯ ¯ ¯ ¯ dt ¯

(6)

with Kc = (16/3π)TE ∆LM AX

Figure 4 Semi-major axis deviation (top), LAN deviation (middle), magnitude of along-track maneuvers (bottom) computed by AOK during the guidance process. Test scenario: constant drag, no other perturbation

The velocity increment to be applied is: µ ∆Vc =

∆ac 2a

¶ V

(7)

where V is the mean velocity. Figure 4 shows an example of guidance process in the ideal case in which the atmospheric density has been fixed at the value 12 · 10−14 kg/m3 (the average of the Harris-Priester model at 700 km altitude) and the atmospheric drag is the only orbit perturbation. At each maneuver cycle of about 1.3 days the evolution of the LAN deviation is evenly parabolic while that of the semi-major is linear. The simulation is performed over 35 days. Figure 5 shows the guidance process in a test scenario in which still the atmospheric density is held constant but Sun-Moon, radiation pressure and tides perturbations are considered. In this case some unevenness is introduced mainly by the effect of the Sun-Moon perturbation action. Atmospheric Drag Estimation By considering Equation (6) it is straightforward to conclude that a correct and fine estimation of the correction maneuver is strictly connected to a correct estimation of da/dt which, as already outlined (Figure 3),

Figure 5 Semi-major axis deviation (top), LAN deviation (middle), magnitude of along-track maneuvers (bottom) computed by AOK during the guidance process. Test scenario: constant drag, Sun-Moon, solar radiation pressure, tides and relativity gravity field perturbations

is mainly determined by the atmospheric drag. As at the altitude of MAIN the atmospheric density has a remarkable variability with time, mainly due to the variability of the solar flux, an atmospheric model cannot give the required accuracy. On the other hand the accelerometers on-board the MAIN satellite are not suitable for the estimation of the atmospheric drag as the typical bias is 4 · 10−2 m/s2 , much larger than typical atmospheric drag acceleration values at 700 km (< 10−5 m/s2 ). The problem has been solved starting by the simple consideration that all the information required for the estimation of the da/dt is available on-board the spacecraft: the knowledge of the actual orbit estimated by the on-board navigation filter with an accuracy of 2 m 1σ (Ref. 9) and the reference orbit (uploaded or computed on-board). The value of da/dt is then found with the following steps: 1. Filter the actual and reference orbit semi-major axis data by taking their values at the ascending node passes. 2. Compute the differences ∆aAN between the actual and reference orbit semi-major axis values. 3. Smoothing of data if required. 4. Find da/dt as the angular coefficient of the line representing the linear fit over the time of a sample of ∆aAN data. Figure 6 shows this procedure applied in the case in which the real orbit is generated by numerical orbit propagation including the aspheric Earth gravity field (GGM01S model) through an expansion in spherical harmonics up to degree and order 90, tides and relativity gravity field perturbations, the Sun and Moon third body forces, the atmospheric drag (atmospheric density model Harris-Priester) and solar radiation pressure. The reference trajectory has been generated using the same gravity field model used to generate the real orbit. Navigation errors are not included. Figure 7 shows the corresponding guidance process in which a sample buffer of 10 ∆aAN ’s is renewed at each ascending node passage for the estimation of da/dt. The buffer is

Figure 6 ∆aAN data fitting process

reset after each maneuver. It can be noticed that the first maneuver is over-estimated because at this point only 5 samples were available and, as no fitting was possible, the da/dt has been estimated using a backup atmospheric model.

Reference Orbit Generation The reference orbit can be generated either on-ground and uploaded on-board the satellite during a ground station pass or directly on-board. The great advantage of the on-ground generation is that a numerical orbit propagation including all the available harmonics of the gravity field model can be used thus achieving the best control performance. The disadvantage is that the reference orbit has to be uploaded periodically (3 days in the case of MAIN) imposing a short time limit on the autonomy of the orbit control. On the other end while an on-board reference orbit generation allows much longer autonomous control periods, the performance is degraded by the fact that only limited gravity field models can be exploited due to the limitation imposed by available on-board CPU resources. In the case that a poor on-board generated orbit is used, the estimation of da/dt turns to be difficult. The short period oscillation of the ∆aAN cannot be catched due to the lack of short period perturbations in the gravity field model. Figure 8 shows the comparison of the ∆aAN curves obtained with a reference orbit generated with a 90x90 gravity field model and one generated with a 5x5 model. In the case of the 5x5 orbit model a short period (for example 10 points) linear fitting is totally unreliable for estimating the atmospheric drag. A long period (50 points) linear fitting would allow an accurate estimation of da/dt, but is not useful at the end of the orbit control because it would require a minimum maneuver cycle period of 50 orbits (i.e. more than 3 days). Nevertheless an efficient estimation of the da/dt is still possible using samples of 10 ∆aAN values but this time applying a proper smoothing process. The smooth filter takes as input the points marked by circles in Figure 8 and flattens them on a line λ keeping unaltered the ratios of their distances from it. Line λ is defined by point (x1 , y) and angular coefficient m where x1 is the x-coordinate of the first point, y is the mean of the y-coordinates of all the points of the sample and m is the estimated da/dt. A guess value of da/dt is used for the very first smoothing process. The output of the smooth filter are the points marked by triangles.

Figure 7 Semi-major axis deviation (top), LAN deviation (middle), magnitude of along-track maneuvers (bottom) computed by AOK during the guidance process (external uploaded reference orbit). Test scenario: drag (Harris-Priester model), Sun-Moon, solar radiation pressure, tides and relativity gravity field perturbations

AUTONOMOUS ORBIT KEEPING The AOK experiment is a minimal target control problem as it has to demonstrate the autonomous fine control of the LAN with respect to a predefined reference orbit. The control accuracy with respect to a predefined phasing grid fixed to a given reference frame on Earth (for example a given longitude) can then be assessed by analyzing the relative movement of the reference orbit and the phasing grid. AOK operations will commence after the end of the PRISMA nominal operations. The test campaign foresees 1 week of initial configuration, calibration and open-loop tests and 4 weeks of close-loop operations. In this section the AOK controller logic and the software development and testing environment are described. AOK Features The AOK control cycle is computed from Equations (6) with da/dt estimated through the process described in the preceding section. The main control logic is the following: 1. An along-track velocity increment ∆V ( Eq. (7) ) is commanded if the time elapsed from last maneuver is larger than a pre-defined time interval ∆tV + , the LAN deviation exceed a pre-defined upper limit (∆LAN > ∆LM AX ) and ∆LAN is growing with time (d∆LAN /dt > 0). Constants ∆LM AX and ∆tV + are uploaded to the spacecraft via Telecommand (TC). 2. After the first along-track maneuver in order to calibrate ∆ac on the actual atmospheric drag value, it is recomputed using Eq. (8): ¯ ¯ ¯ da ¯ Tac ∆ac ∆acal = ¯¯ ¯¯ + dt ac 2 2

(8)

Figure 8 Estimation of da/dt for a reference orbit generated with a 5x5 gravity field model

where: (da/dt)ac is estimated by the long period fitting (with the number of points imposed by TC) if possible or by the ordinary short period fitting, and ∆acal is the calibrated value of ∆ac . Tac is the actual cycle period i.e. the time elapsed from the actuation of last maneuver. 3. An anti-along-track velocity increment is commanded if this possibility is enabled by TC, the LAN deviation exceeds a pre-defined lower limit (∆LAN < −n·∆LM AX where n is a real constant), d∆LAN /dt < 0 and the time elapsed from last maneuver is larger than a pre-defined time interval ∆tV − . Constants n and ∆tV − are uploaded to the spacecraft via TC. In this case: ∆acal

¯ ¯ ¯ da ¯ = ¯¯ ¯¯ dt

ac

Tc B

(9)

where: Tc is computed with Eq. (6) B is an empirical correction factor imposed by TC. Figure 7, provides an example of the AOK guidance process in case of a reference orbit is uploaded from ground. ∆LM AX is 10 m, the anti-along-track maneuvers are not enabled, ∆tV + is zero. As the reference orbit is the best available and neither navigation nor sensors and actuators errors are introduced in the simulation, ∆LAN has a smooth parabolic evolution. The mean of ∆LAN is -5 m and σ (standard deviation) is 12.3 m: the control accuracy of the osculating ascending node is then compliant with the main requirement (10 m 1σ). In contrast Figure 9 gives an example of AOK guidance with an on-board propagated orbit (gravitational field model GGM01S 5x5). ∆LM AX is 50 m and ∆tV + is 5 orbital periods (about 8 hours). The enlarged control window and waiting times are chosen to compensate the presence of oscillations in the evolution of ∆LAN which are due to the incompleteness of the reference orbit. The enlargement of this two parameters allows AOK to correctly estimate da/dt. A parabolic but noisy evolution of ∆LAN can still

Figure 9 Semi-major axis deviation (top), LAN deviation (middle), magnitude of along-track maneuvers (bottom) computed by AOK during the guidance process: on-board generated reference orbit, ∆LMAX = 50 m, ∆tV+ = 0

Figure 10 Semi-major axis deviation (top), LAN deviation (middle), magnitude of along-track maneuvers (bottom) computed by AOK during the guidance process: on-board generated reference orbit, ∆LMAX = +50 m ÷ -100 m, ∆tV+ = 8 h, ∆tV− = 0

be observed. This time the mean of ∆LAN is -50 m and σ is 53 m: the control accuracy is not compliant with the main requirement. Figure 10 shows the same test case as Figure 9 but with the enabling of anti-along-track maneuvers: the control window is this time bounded between -100 m and +50 m ( ∆LM AX = 50 m, n = 2). In this case the performance of AOK is better (the mean of of ∆LAN is -22 m and σ is 45 m) but at the cost of a larger number of maneuvers which impose a negative ∆V most of the times. Software Development and Testing Environment Figure 11 shows the simplified architecture of the PRISMA on-board software (OBS) architecture for the MAIN satellite and the simulation environment (see also Ref. 18). The PRISMA OBS has a layered structure with a Basic Software (BSW) level and an Application Software (ASW) level communicating with each other through dedicated message queues. While the BSW includes basic applications, device drivers and I/Outilities, the ASW encapsulates all top-level applications like spacecraft navigation, control, telecommand and telemetry. The GPS-based navigation system is split into three modules located in different OBS levels and running at different sample rates. The GPS interface (GIF) is part of the BSW, runs at 1 s sample time and is directly fed with GPS messages issued by the Phoenix GPS receivers emulator (PEM). GIF handles GPS raw data formats and ephemerides, and performs data sampling as well as coarse editing prior to the GPS-based orbit determination. The GPS-based Orbit Determination (GOD) and GPS-based Orbit Prediction (GOP) are embedded in the ASW layer as part of the ORB core (30 s sample time) and the GNC core (1 s sample time), respectively. GOD implements an extended Kalman filter to process GRAPHIC observables as well as single difference carrier phase measurements from MAIN satellite. Attitude data from both spacecraft are applied to correct for the GPS receivers antenna offset with respect to the spacecraft center of mass. Furthermore, a history of maneuver data is provided to GOD and taken into account in the orbit determination task. GOD performs a numerical orbit propagation which is invoked after the measurement update and provides orbit coefficients for interpolation to GOP for both spacecraft. The GOP module interpolates the orbit coefficients provided by GOD and finally supplies the various GNC core functions as well as the PRISMA payload with continuous position and velocity data. Due to the different data rates of the GPS-based navigation modules, orbit maneuver data have to be taken into account in both GOD and GOP (cf. Figure 11). In particular at each GNC step, the GOP task accounts for maneuvers which have not been considered by GOD in the last orbit determination/prediction process. The AOK and Autonomous Formation Control (AFC) modules are located in the ORB and GNC cores respectively. Finally the orbit propagator (PROP) generates the MAIN simulated orbit. For PRISMA, a model-based design method is used for the complete on-board application software. While traditional software development is performed on a programming level and focused on software engineering, the model-based approach emphasizes systems engineering methods.

Figure 11 Simplified scheme of the simulation software architecture

This raises the abstraction level for the development in order to allow an efficient handling of ever more complex systems. Benefits of a model-based over a traditional development approach are expected in the field of robustness, predictability, scalability, efficiency and faster development schedules. The software is implemented in Matlab/Simulink as well as C/C++. In particular, C++ low level routines, which are interfaced with Matlab/Simulink high level structures, represents the computational layer of the software system including for example numerical integrations, data filtering, etc. Matlab/Simulink provides instead the communication layer, including for example input/output interfaces, complex logics, time synchronization. The prototyping of the software as well as the related analysis and simulations are performed first on a host standard laptop PC in a Matlab/Simulink environment. All the functions can then be auto-coded using the Matlab tool Real Time Workshop and executed under the operating system Real Time Executive for Multiprocessor Systems (RTEMS) on the target LEON3 processor (the PRISMA satellites on board computer). The full consistency of the flight software generation chain can thus been verified by comparing relevant outputs generated by the software running on the host and target computer. One of the greatest advantages of this straightforward and cheap approach is that it allows the test of the flight software performance directly on the on-board spacecraft computer. Table 3 AOK execution times and CPU loads on LEON3 board AOK Mode

Execution Time (ms)

CPU Load (%)

63 596 49623

0.2 2 165

External uploaded RO On-board generated RO (multi-step) On-board generated RO (one-step)

The analysis of the software execution times has special interest because it allows to verify if the different software modules comply with the maximal allocated CPU loads. Analysis of heap and stack allocation and software profiling can also be made for the detection of eventual software execution bottlenecks. Table 3 collates the execution times and CPU loads of AOK running on the LEON3 board in three modes: with an external uploaded RO, a RO generated on-board by a propagation performed in one call or in several calls. Only the one-step propagation is not compliant with the CPU load constraints. External uploaded and multi-step on-board propagated RO modes are fully compliant with the CPU load constraints because AOK is located in the ORB core that is a low priority task. The possibility of using AOK in both external and internal propagated RO modes, enhances the robustness and the autonomy of the control as AOK can automatically switch to on-board generated RO in case of not validity of the external uploaded RO. Numerical Simulation Figure 12 shows how the guidance process of AOK in the test scenario of Figure 7 changes if the navigation errors are included and sensors and actuators are modeled by applying realistic errors to measurements and control actions respectively. Table 4 collects the navigation accuracy and the parameters used for the accuracy modeling of sensors and actuators. The propulsion system is characterized by a Minimum Impulse Value (MIV) of 0.7 mm/s and a Minimum Impulse Bit (MIB) of 0.07 mm/s. Consequently, the thrusters can only realize ∆V ’s which are larger than MIV and integer multiples of MIB. Furthermore, the performance error of the thruster system is quantified by Equation (10) ¯ ¯ ¯ ∆Vr − ∆Vp ¯ ¯ ¯ · 100 ²=¯ ¯ ∆Vp

(10)

where ∆Vp is the planned velocity increment computed by the AOK controller and ∆Vr is the actual velocity increment realized by the propulsion system. Finally the attitude control error, which causes deviations from the desired thrust direction, is treated as Gaussian noise with zero bias and 0.2◦ standard deviation in three axes. The following conclusions follow from an analysis of the simulation results shown in Figure 12:

Figure 12 Semi-major axis deviation (top), LAN deviation (middle), magnitude of along-track maneuvers (bottom) computed by AOK during the guidance process

1. The inclusion of navigation errors causes an oscillating evolution of ∆a and consequently a different evaluation of the LAN deviation with respect to the case in which no errors are included. It follows that the maneuver planning is different in the two cases. 2. The inclusion of navigation errors decreases the estimation accuracy of da/dt. It follows that, as the efficiency of the commanded maneuvers is reduced, a larger number of thrust pulses is needed in the guidance process. 3. The inclusion of sensors and actuators errors further reduces the efficiency of the commanded maneuvers. 4. In the case in which all the errors are included the mean of ∆LAN is -5 m and σ is 17 m. Table 4 Navigation, sensors and actuators accuracy modeling Absolute Navigation Accuracy

Value

Mean σ

0 (m) 2 (m)

Attitude Control Accuracy

Value

Mean σ

0 (◦ ) 0.2 (◦ )

Propulsion System Accuracy MIV MIB ²

Value −4

7 · 10 (m/s) 7 · 10−5 (m/s) 10

CONCLUSION The AOK experiment is the secondary objective of the DLR contributions to the PRISMA formation flying mission. It has to demonstrate autonomous orbit control of the longitude of the ascending node of a single spacecraft with an accuracy of 10 m (1σ). The software implements an analytical feedback control algorithm using along-track and anti-along-track velocity increments. In the realization of the AOK controller the problem of a correct estimation of the atmospheric drag has been solved. The software applies a filtering, linear fitting and smoothing of the difference between the reference and real orbits’ semi-major axis values at the ascending node passes. Advantages and disadvantages of using an on-ground or on-board generated reference orbit have been analyzed in terms of performance and operativeness. The results of the simulation performed are encouraging in showing that the fulfillment of the strict control accuracy requirement is feasible. Execution time tests on a representative PRISMA on-board computer indicate that at the present state the software is fully compliant with the CPU load requirements.

REFERENCES [1] J. R. Wertz, “Autonomous Navigation and Autonomous Orbit Control in Planetary Orbits as a Means of Reducing Operations Cost,” 5th International Symposium on Reducing the Cost of Spacecraft Ground Systems and Operations, Pasadena, CA USA, July 2003. [2] S. Kudryavtsev and E. Gill, “Autonomous Orbit Control of Remote Sensing Spacecraft using Spaceborne GPS Receivers,” 14 December 2000. [3] S. D’Amico, C. Arbinger, M. Kirschner, and S. Campagnola, “Generation of an Optimum Target Trajectory for the TerraSAR-X Repeat Observation Satellite,” 18th International Symposium on Space Flight Dynamics, Munich, Germany, 11-15 October 2004. [4] C. Arbinger, S. D’Amico, and M. Eineder, “Precise Ground-In-the-Loop Orbit Control for Low Earth Observation Satellites,” 18th International Symposium on Space Flight Dynamics, Munich, Germany, 11-15 October 2004. [5] M. Rosengren, “The Orbit Control of ERS-1 and ERS-2 for a very accurate Tandem Configuration,” RBCM, Vol. XXI Special Issue, 1999, pp. 72–78. [6] M. Rosengren, “Improved Technique for Passive Eccentricity Control,” AAS 89-155, AAS Symposium on Mission Design and Orbital Mechanics, 1989. [7] S. Persson, P. Bodin, E. Gill, J. Har, and J. J¨orgensen, “Prisma An Autonomous Formation Flying Mission,” Sardinia, Italy, ESA Small Satellite Systems and Services Symposium (4S), 25-29 September 2006. [8] E. Gill, S. D’Amico, and O. Montenbruck, “Autonomous Formation Flying for the PRISMA Mission,” Journal of Spacecraft and Rockets, Vol. 44, May 2007, pp. 671–681. [9] S. D’Amico, E. Gill, M. Garcia, and O. Montenbruck, “GPS-Based Real-Time Navigation for the PRISMA Formation Flying Mission,” Noordwijk, 3rd ESA Workshop on Satellite Navigation User Equipment Technologies - NAVITEC, 11-13 December 2006. [10] A. Lamy, M.-C. Charmeau, D. Laurichesse, M. Grondin, and R. Bertrand, “Experiment of autonomous Orbit Control on the Demeter Satellite,” 2004. [11] A. Lamy, M.-C. Charmeau, and D. Laurichesse, “Experiment of autonomous Orbit Control on the Demeter Satellite: In-flight Results and Perspectives,” AIAA Astrodynamics Specialist Conference, Keystone, Colorado, USA, 21-24 August 2006. [12] xPC Target 3 - Users Guide. The MathWorks, March 2008. [13] J. Gaisler, The LEON Processor User’s Manual V2.3.7. Gaisler Research, 2001. [14] RTEMS C Users Guide. Edition 4.6.5, for RTEMS 4.6.5. On-Line Applications Research Corporation (OAR) http://www.rtems.com, 30 August 2003. [15] J. Liu, “Satellite Motion about an Oblate Earth,” AIAA Journal, Vol. 12, November 1974, pp. 1511–1516. [16] O. Montenbruck and E. Gill, Satellite Orbits - Model, Methods and Applications. Heidelberg: Springer Verlag, 2000. [17] P. Micheau, Spaceflight Dynamics, Vol. I, ch. Orbit Control Techniques for LEO Satellites. Toulouse, France: Cepadues-Editions, 1995. [18] S. D’Amico, S. D. Florio, and T. Yamamoto, “Offline and Hardware-in-the-loop Validation of the GPS-based Real-Time Navigation System for the PRISMA Formation Flying Mission,” 3rd International Symposium on Formation Flying, Missions and Technology, ESA/ESTEC, Noordwijk, The Netherlands, 23-25 April 2008.