Relative State Transition Matrix Using Geometric Approach for Onboard Implementation

: In the present study relative state transition matrix was obtained. It matches the relative motion of 2 satellites while including the oblate perturbation. The used formulation applies the geometric approach and is in Cartesian frame. The relative state transition matrix uses absolute state transition matrix of individual satellites. Thereon, simplification in computing methods and on-board implementation at controls are explored in a leader/follower coordination method. Numerical experiments illustrate the accuracy when the baseline separation is equal to 2 km and eccentricity is 0.005 and 0.05 for all the inclinations.


INTRODUCTION
Formation flying is about harnessing and controlling the relative orbital dynamics.When using multiple observations obtained from a formation of the low Earth orbiting satellites, the combined data is enhanced.The Clohessy-Wiltshire (CW) equation that describes a relative motion is valid for a small relative distance between the satellites compared to the radius and it assumes that the reference orbit is circular besides using a spherically uniform gravitational field at the centre (Alfriend et al. 2010).In practice this is not true for a low Earth orbiting satellite.Maintaining a formation with respect to the circular reference orbit requires prohibitive fuel consumption.It needs frequent orbit corrections to reduce the deviations caused by perturbations.The first significant perturbation is due to the effect of oblate Earth denoted as J 2 .When the reference orbit is initially nominally circular, a small eccentricity usually occurs due to the J 2 effect, which also causes a regression in the orbital plane.It may be noted that an advanced scheduling of the satellite payload data collection on the Earth is usually based on a ground system.This mapping system uses an orbit model that includes mean oblate effect.Thus, in general, it is more important to include the J 2 effect in the on-board reference orbit of the satellite.
State transition matrix (STM) plays an important role in the orbit control.It occurs in the system dynamics equation that is used to control the satellite.The STM that is implemented on-board a satellite needs to be simpler, computationally efficient, and yet scaled close to match the accompanying reference on-board orbit model.Here, the reference orbit on Relative State Transition Matrix Using Geometric Approach for Onboard Implementation each of the satellites is a first-order motion that includes the J 2 effect (Roy 1982).In this paper, it was obtained an approximate computationally simpler relative STM in Cartesian frame to match this relative motion.
There are various formulations of relative STM and each assumes certain orbital characteristics (Alfriend et al. 2010).Here a rectangular Cartesian frame is used since it helps in the direct use of GPS navigation systems data.Incrementing the maneuver velocity is easier to understand with respect to the orbit frame of reference.Cartesian frame has been found to perform satisfactorily under propagation for formation types that are not largely separated along-track (de Brujin et al. 2001).
The first STM is that of the circular CW equation.The assumptions in the CW model, as mentioned earlier, are the limitations in its use.The STM is derived in Cartesian frame from the analytical solution of the CW differential equation (Alfriend et al. 2010).When the eccentricity of the reference orbit is small, a STM can be derived in Cartesian frame, as carried out by Melton (2000).This involves a lengthy expression as a function of time and is a series with terms which are powers of the eccentricity (Melton 2000).When the orbit is more elliptical then the Keplerian one, it is transformed following the approach of Tschauner and Hempel (1965) to yield a closed form relative motion.This STM given in Yamanaka and Ankersen (2002) is an improvement made to Carter (1998).Here the independent variable is the true anomaly.The time-explicit STM is easier to adopt in the control system.
The method of Gim and Alfriend (2003) has solved the problem of STM, having matched the relative motion about an elliptical orbit with J 2 effect.Though the time-explicit expression is analytical and accurate, it is laborious for on-board implementation.The model includes short and long periodic effects and is in terms of non-singular elements.It is important to mention that there is no need to include the long-period effects under J 2 perturbation.In Hamel and de Lafontaine (2007), this is explored, and an alternative STM was presented, being valid for all cases when the reference orbit is not circular.It is based on mean orbit element differences between 2 satellites.When the eccentricity is small, then a lesser complex version is given in Alfriend and Yan (2005).The STM formulated in Yan et al. (2004) is based on the unit-sphere method (Vadali 2002).The approximate STM in Tsuda (2011) is in Cartesian frame and involves a series for on-board implementation.
Here the approach is geometric, that is, it applies a geometric transformation between the relative states of the deputy with respect to the chief and the orbital frame of the chief satellite.This was adopted to describe the relative motion in Balaji and Tatnall (2003) which is similar to the unit-sphere approach (Vadali 2002) and appeared at the same time.The relative motion described in Vadali (2002) as well as in Balaji and Tatnall (2003) includes the J 2 effect, is applicable for eccentric orbits, and is valid for larger relative distances, using Keplerian elements.The present geometric approach is in Cartesian frame.The relative STM obtained here is approximate and explores optimized computation.Navigation plan suggested here is suited for on-board real-time usage of STM in leader/follower coordination (Alfriend et al. 2010).The STM is applicable for both osculating and mean elements of the 2 satellites.It also provides a scope of including secular effects due to higher-order perturbation.Finally, numerical experiments bring out the accuracy of state updated by the STM in terms of varying eccentricity and inclination.

REFERENCE ORBIT AND STATE TRANSITION MATRIX
The orbital motion of the deputy and chief satellite with J 2 effect is discussed next.The absolute STM of the 2 satellites is then derived.This is required in each of the satellites in controlling their reference orbits.
The equation of motion of an individual satellite is (Roy 1982): where R ¨ denotes the second derivative with respect to time of R = (X,Y,Z).The disturbing potential is: (1) (2) (3) with where: μ is the Earth's gravitational parameter, equal to 398,600.64 km 3 /s 2 ; J 2 = 0.0010826; δ is the instantaneous declination; R e is the Earth's radius, which is 6,378.458km; R is the magnitude of the satellite position vector, R.
The position vector of the deputy and the chief satellites is denoted by R d and R c , respectively, in the Earth Centred Inertial Frame (ECIF).The potential is axi-symmetric about the Z-axis and is independent of the azimuth angle ψ.
F has 2 types of expressions (Roy 1982): for close formations to avoid any collision.On the other hand, having the reference to be a mean motion fuel is saved by having impulse maneuver used once or twice appropriately in the orbit besides maximizing data collection duration (Gill et al. 2007).
) be in orbital state in the ECIF.At the epoch (t 0 ), the state is X j (t 0 ), where j = 1 denotes the state of the chief and j = 2 that of the follower.
The absolute STM formulation follows the method of Markley (1986).Equation 9describes the potential and is in Cartesian frame when using periodic and secular terms: where: F 1 contains only the secular component while F 2 includes the short-period perturbation.
where f is the true anomaly.The Lagrange's planetary equation of motion (Roy 1982) is: where: where: a, e, and i are, respectively, the semi-major axis, eccentricity, and inclination; ω, Ω, and M are, respectively, argument of perigee, longitude of the ascending node and mean anomaly.
There are 2 kinds of reference orbits for on-board implementation on both satellites.It is either (a) osculating and obtained by substituting F 2 in Eq. 6 or (b) the mean motion obtained by substituting F 1 in Eq. 6.The latter is the more familiar first-order secular perturbation with a, e, and i remaining constant, as it is at the start epoch.When mapping the orbit elements to the mean ones, only 3 orbit elements are found to exhibit secular growth due to the J 2 effect and are where: R j = (X j 2 + Y j 2 + Z j 2 ) 1/2 .On the other hand, the potential for secular effects alone is: The acceleration derived from the potential for the case (a) with secular and short periodic motion is presented as in Chairadia et al. (2012): On the other hand, for (b), when only secular effects are considered, the accelerations can be seen in Ramachandran (2015) as: It can be noted that the inclination i is invariant in the motion described in Eq. 8.These accelerations are then used to derive the partial derivative and obtain the absolute STM matching the 2 types of reference motions: (a) or (b).Details are given in Markley (1986), as well as in Chairadia et al. (2012).

It is more meaningful to use reference orbit with periodic terms
Relative State Transition Matrix Using Geometric Approach for Onboard Implementation in the chief satellite's local-vertical-local-horizontal coordinate frame (LVLH).The relative position vector r, with respect to the chief, has 3 components: x, y, and z, as illustrated in Fig. 1.Here, (x) is the component in the chief 's radial direction, (z) is in the direction of the chief instantaneous angular momentum, and (y) is along the direction that completes the triad (Alfriend et al. 2010).The relative distance is: Here [P] is a (6 × 12) matrix.In case the mean elements of the state Š are used, then another transformation needs to be carried out as in Gim and Alfriend (2003), in addition to that of using Cartesian transformation.for the chief and deputy satellites, respectively.The instantaneous angular momentum vector of the chief is H c = (R c × V c ), where (x) is the usual vector cross-product.In Eq. ( 13.2), the second row is the normalized cross-product of the first and third rows.
Taking the derivative of Eq. 13.1, the relative velocity is: The transformation matrix D is reversible and has been sufficiently discussed by Gim and Alfriend (2003).
The following discussion uses osculating element.The formulation is analogous while using mean elements after applying the transformation as in Eq. 16.The relative STM Φ r (t,t 0 ) follows: with The matrix P is known in Eq. 15.The second term on the right-hand side of Eq. 18 is the matrix Q and it uses the absolute STM derived earlier for both satellites.If the input is the reference state S(t 0 ), then the update is obtained by using Using Φ a to individually denote the STM of the chief and deputy satellites, we have The matrices are the absolute STM obtained either (a) using periodic accelerations or (b) admitting only secular accelerations as mentioned in Eqs.11 and 12, respectively.In case the relative STM for mean elements is needed, then the transformation for the osculating to mean, as in Eq. 18, is required.
In order to derive the elements of C oi 1 in P, a general relation is used x rel = [P ] S x rel = [P ] S x rel = [P ] S Combining to get Substituting the 3 rows of C oi in place of (u/|u|) in the above equation, it is obtained x x where Ramachandran MP (21) Again the second row is computed as a cross product of the first and third rows instead of the algebraic expression.

COmpuTATiOnAl REDuCTiOn
The computation reduction aspect is discussed next.The absolute STM are anyway needed on-board for individual orbit control.The additions for relative control are the matrices [P] and [T] and subsequently the matrix multiplication.It is now examined how this computation of the relative STM in Eq. 18 can be decomposed for computation efficiency.It may be noted that there are 4 matrices.
The chief or target reference orbit is available in the deputy or follower satellite due to inter-satellite communication (Gill et al. 2007).The thruster on the deputy then controls the orbit to maintain the formation in a leader/follower coordination (Alfriend et al. 2010).The chief reference orbit is undisturbed during this operation.Using the states at t 0 and t 1 , the elements of the matrix in Eq. 25 can be computed, ahead of that time.This management ensures that the elements of the first 3 rows of the matrix can be obtained beforehand.The orbit of the deputy or main satellite and the matrix [λ] alone varies due to orbit control thrusters and needs to be updated.This computational reduction in real time is the advantage of the present approach.On the other hand, in the case of unit-sphere approach (Vadali 2002), the transformation matrix in Eq. 15 is constructed using both chief and deputy states.Hence, all the elements of the relative STM need to be computed in real time and are task-intensive.The STM in Gim and Alfriend (2003) is more difficult for implementation.
While using mean motion, the secular accelerations given in Eq. 11 are used to construct mean relative STM.This, however, calls for conversion of the measurement transformed to first-order mean elements using Eq.16 and that further to Cartesian frame.The orbit model in Eq. 8 can be enhanced to accommodate secular effects of other perturbation for near circular orbits.Towards this, the orbit model suggested in Ramachandran (2015) can be adopted.It means that the motion is more realistic to account for secular variation which is along the track between the satellites in formation, and the STM is useful to match the realistic motion of the satellite.

NUMERICAL ANALYSIS
The experiment carried out here intended to estimate the accuracy of the STM given in Eq. 17.The mean motion is selected as input for each equation.From Eq. 8, the initial states at times t 0 and t 1 are achieved, and the matrix [Q] is obtained.Using the known relative position at time t 1 as input and the approximate STM in Eq. 17, the update of the relative state at time t 2 is obtained.This appears on the left-hand side of Eq. 17.The time step selected is 1 s.Subsequently, it is obtained the difference between this output ∆x rel (t 2 ) and the actual relative state at that instant computed using Eqs.13 and 14.The error in position and velocity in the LVLH frame is related to the reference orbit.This error is gauged against variations in inclination and eccentricity.It is known that the J 2 effect brings an out-of-plane (z) effect which is absent in CW equations.In the first experiment, the orbital eccentricity is 0.05, and the semi-major axis a is 7,000 km.The deputy satellite is at a distance of 2 km along the track.The plot in Fig. 2 shows the error while varying inclination from 0 up to 90 degrees when retaining the eccentricity as 0.05.The numerical accuracy of the methods from Gim and Alfriend (2003) and Yan et al. (2004) shall be no different as they are mathematically similar and differ in formulation.Also, the comparison of propagation error will be limited by the frame of reference.A case of higher eccentricity when e = 0.1 is carried out.The error in the position and velocity was observed to be greater than that shown in Fig. 2, though the variation with inclination is similar.The errors in position and velocity are −0.68 m and −0.41 m/s, respectively.This increase in the error is as expected for higher values of e and due to the approximation made in Q of Eq. 19.2.
It is seen that x and x 1 (see Eqs. 13 and 14) have more perceivable error than (y, y 1 , z, z 1 ).The position error is in the top plot while the velocity error is that at the bottom.The (y, z) values agree to 4 decimal places while (y 1 , z 1 ) ones agree beyond 5 decimals.This implies that the update error along these axes is negligible.
To understand the effect of eccentricity, the experiment is repeated by changing the eccentricity to 0.005.The units now are in cm and cm/s -yet again the error is large in (x, x 1 ).The error in the other axes agrees in 2 decimal places in position and 3 in velocity.When the distance of the satellites is 10 km, then the error gets proportionally scaled to that in Fig. 3.This is due to the approximation in the absolute STM.
The next test is carried out by varying the eccentricity from 0.0006 to 0.06.The inclination is fixed and equal to 20 degrees.In Fig. 4, the error along (x, x 1 ) is greater again and increases with eccentricity.It can be seen in Markley (1986) that the error in absolute STM depends upon the position in the non-circular orbit.In this example, the argument of perigee is retained while the eccentricity varies.
As a summary, it is clear that the error, using Eqs.13 and 14, is the least along the track and across the track (y, z) -these are very important in the control of the satellite.The radial error along (x, x 1 ) was found more agreeable for smaller eccentricity.

CONCLUSION
Based on the geometric approach, a relative STM was obtained, and it includes the oblate Earth effect.It is simpler for on-board implementation as it is decomposable in terms of the orbital states of the 2 satellites in formation.The navigation plan, as suggested in this paper, reduces the real-time computation when applying the proposed relative STM for orbit control purposes.A numerical experiment was carried out to measure the accuracy of the proposed method for all inclinations.It shows that the approach is more accurate when the distance between the satellites is below 2 km and when the eccentricity is 0.005.For larger baselines and eccentricity, the accuracy of the relative STM is seen limited due to the first-order approximation in the absolute STM.
C oi (3 × 3) relates the LVLH frame of the satellite to the ECIF.The position and velocity vectors are, respectively, R c = (X j ,Y j ,Z j ) and V c = (j = 1 denotes the chief satellite while R d and V d for j = 2 are used to denote the deputy satellite.The magnitudes are |R c | and |R d |

Figure 1 .
Figure 1.Relative motion of deputy satellite in the frame of the chief.

Figure 4 .
Figure 4. Error variation with respect to eccentricity when inclination is 20 degrees.