# Approximate State Transition Matrix and Secular Orbit Model

Authors: M. P. Ramachandran
Journal: International Journal of Aerospace Engineering
 (2015)

---

## Abstract

The state transition matrix (STM) is a part of the onboard orbit determination system. It is used to control the satellite’s orbital motion to a predefined reference orbit. Firstly in this paper a simple orbit model that captures the secular behavior of the orbital motion in the presence of all perturbation forces is derived. Next, an approximate STM to match the secular effects in the orbit due to oblate earth effect and later in the presence of all perturbation forces is derived. Numerical experiments are provided for illustration.


---

## Body

## 1. Introduction

Autonomous orbit control in satellites is possible with the present onboard technological advancements. The Global Positioning System receiver solution gives the satellite position measurement in Cartesian frame. State propagating equations along with the measurement equations in the linear filter then estimate the orbit. State transition matrix (STM) is used in the state update equations. A reference orbit model is available onboard. Using the receiver orbit solution the absolute orbit control system then ensures the satellite motion to this reference orbit model in the earth centered fixed Cartesian reference frame. This control enables the satellite to achieve the required orientation too. In orbit determination system, STM of two-body dynamics as suggested in [1] is usually used. Yet it will be always desirable to match the complete dynamics especially to improve the accuracy and scalability of the navigation system [2].The orbital motion of the satellite is made up of secular or mean motion along with short and long periodic motions [3, page 571]. When we include the complete dynamics as reference orbit we have to use continuous control. This requires more fuel. Continuous maneuver can also disturb the payload functioning. On the other hand mean motion (without periodic motions) as a reference orbit is more suited for orbit control by impulse thrusting. This is adopted in formation flying [4, Chap 10] besides that the mean motion is used to derive the initial conditions. It is noted that orbit control is usually executed as a function of time [5] instead of true anomaly. In the control system, the state measurements in Cartesian frame are usually updated in time space. Subsequently the STM derived here then updates these states.STM henceforth shall mean absolute STM unless mentioned. It may be noted that Vallado [3, page 748] has discussed the STM for two-body orbital motion. In [6] a STM including the oblate earth effects using equinoctial mean elements and then applying interpolation is obtained. The present note brings out a STM that is in Cartesian frame as an alternative to [6] and considers only secular effects. We note in the literature that the STM that is in Cartesian frame is derived in [7] and it includes secular and periodic effects. Here the periodic effects are neglected. Further the STM derived here is extendable to accommodate secular along track effects in the presence of all perturbations. This is simpler than the expansion based method of deriving the STM as in [8].It is important to note that secular forces due to oblate effect are considered in relative motion as in formation flying which is based on geometric approach [9]. The absolute transition matrix derived here can further be used to derive relative transition matrix as in [10]. This work is beyond the scope of this paper.


## 2. Secular Acceleration

Consider the equation of motion(1)r¨=∇U+ap,where r¨ denotes the second derivative with respect to time of r=(x,y,z), the position vector in the inertial frame. The disturbing potential [11] is(2)Ur,ψ,z=μr-μrRr2J2123sin2⁡δ-1,where J2=0.0010826, δ is the instantaneous declination, R is the radius of the earth, μ is the gravitational constant, and r is the magnitude of the position vector, r. The vector ap=(apx,apy,apz) represents other perturbation forces due to the inhomogeneous mass distribution of the earth, third body forces due to sun and moon besides solar radiation pressure and atmosphere drag forces. The potential is axisymmetric about the z-axis and is independent of azimuth angle ψ. The Lagrange’s planetary equation of motion is invoked and the following relations are deduced. The Keplerian elements are averaged over an orbit. The first-order secular motion that neglects periodic effects is described by(3a)a=a0,e=e0,i=i0,where a,e,i are, respectively, the semimajor axis, eccentricity, and inclination are invariant over the duration of interest:(3b)w=w0+32J2Rp2nˇ2-52sin2⁡iΔt,Ω=Ω0-32J2Rp2nˇcos⁡⁡(i)Δt,nˇ=n01+32J2Rp21-32sin2⁡i1-e21/2,n0=μa3,M=M0+nˇt.And w,Ω,M are argument of perigee, the longitude of ascending node, and mean anomaly, respectively. The equation of the centre enables getting the true anomaly (f):(4)f=M+2e-14e3sin⁡(M)+54e2sin⁡(2M)+1312e3sin⁡(3M)+o(e4).In (3a) and (3b) we note that(5)p=a(1-e2).Equation (4) is used when the eccentricity is not large. The longitude of the ascending node Ω0 varies linearly with incremental time, Δt. Orbit models in the satellite for control purposes need to have the cross-track motion that is predominant due to J2. The argument of perigee (w) along with the true anomaly (f) gives the argument of latitude (θ) which is(6)θ2=w2+f2.Here, the subscript 2 has been added to denote the J2 model. When all forces are included and solution of (1) is obtained, the instantaneous argument of latitude is denoted by θa. Here, in this paper, a proposal is made to add a polynomial function to the argument of latitude θ2, in (6) over every orbit. This is a mean variation of the differential argument of latitude and could be, say, a quadratic or cubic power of time and is denoted as θf [3, page 570, 652]. Here, in this note, θf is a least squares fit over one orbital period and it accommodates the secular difference. This is defined as(7)θp≝θ2+θf.The residue between θa and θp is periodic, which is incidentally not required for control. This correction enables (3a) and (3b) along with (7) to match the secular effects when all perturbation is present to a reasonable accuracy specially along track. Next a STM that matches the orbit model in (3a) and (3b) and later the secular effect in the presence of all forces is derived.The following equations are used to transform from the orbital frame[r,t,n] to the Cartesian frame [x,y,z]:(8)xˇyˇzˇ=Arˇtˇnˇ(9)a1,1=cos⁡⁡(θ)cos⁡⁡(Ω)-cos⁡⁡isin⁡Ωsin⁡(θ)a1,2=-sin⁡θcos⁡⁡Ω-cos⁡⁡isin⁡Ωcos⁡⁡(θ)a1,3=sin⁡(i)sin⁡(θ)a2,1=cos⁡⁡(θ)sin⁡(Ω)+cos⁡⁡(i)cos⁡⁡Ωsin⁡(θ)a2,2=-sin⁡θsin⁡Ω+cos⁡⁡(i)cos⁡⁡Ωcos⁡⁡(θ)a2,3=-sin⁡icos⁡⁡(Ω)a3,1=sin⁡(θ)sin⁡(i)a3,2=sin⁡(i)cos⁡⁡(θ)a3,3=cos⁡⁡i,where unit vectors rˇ,tˇ,nˇ, respectively, are in the radial, tangential, that is, along the direction of motion (along-track) and normal to the orbital plane (see Figure 1).Figure 1
Orbit frame from position, tangential, and normal vectors.By usingsin⁡(δ)=rsin⁡(θ)sin⁡(i) in (2), the potential due to oblate earth effect becomes(10)U2=-μrRr2J212{3sin2⁡isin2⁡(θ)-1}.The acceleration is then obtained using the relations in (8) and is(11)∇U=-μr2rˇ-μR2J212-15z2r6+3r4(rˇ)-μR2J2126zr5(zˇ).Equation (11) when substituted into (8) finally gives the accelerations in the (x,y,z) frame:(12)ax=-μx(r3)1+3/2J2R/r21-5z2/r2,ay=yxax,az=-μz(r3)1+3/2J2R/r23-5z2/r2.Partial differentiation with respect to the state variables (x,y,z) in (12) yields the STM that includes both secular and periodic components (see [7]).Here, the termcos⁡⁡(2θ) in (10) is periodic and neglected; the remaining term that is secular is retained. In (2) we have(13)U(r,δ)=μr-μrRr2J21232sin2⁡i-1.Note that the potential in (13) is independent of (z). The net acceleration then is obtained as in (11) by denoting α=-μR2J2(3/2)(1-3/2sin2⁡(i)) in Cartesian frame and is given by(14)ax=αxr5-μxr3,ay=αyr5-μyr3,az=αzr5-μzr3.Then (15)∂ax∂x=-3αr5+15αx2r7+3μx2r5-μr3,∂ax∂y=15αxyr7+3μxyr5,∂ax∂z=15αxzr7+3μxzr5,∂ay∂x=∂ax∂y,∂ay∂y=-3αr5+15αx2r7+3μy2r5-μr3,∂ay∂z=15αyzr7+3μyzr5,∂az∂x=∂ax∂z,∂az∂y=∂ay∂z,∂az∂x=-3αr5+15αz2r7+3μz2r5-μr3.It may be noted that when θ is substituted into (8) by θ2 in (3a) and (3b) or θp in (7), the accelerations derived in (15) are still valid. This implies that the STM (to be derived in the section) with the accelerations derivatives in (15) are valid for (1) without considering periodic effects.


## 3. Approximate STM

Next use the total acceleration in (14) and derive the approximate STM following Markley [12]. The STM is then obtained approximately based on Taylor expansion(16)Φt,t0=ΦrrφrvΦvrφvv.With a knowledge of initial states as in (3a) and (3b) at t0, the matrix Φ can be used to obtain state at the subsequent instant “t,” using(17)Xt=Φt,t0Xt0,where X(t) is differential of the states (x,y,z,x′,y′,z′) at t0. Discarding higher order terms in (16) we have(18)Φrr=I+2G0+GΔt26,Φrv=IΔt+2G0+GΔt312,Φvr=G0+GΔt2,Φrv=I+G0+2GΔt26.The gradient matrix is(19)G=∂ax∂x∂ax∂y∂ax∂z∂ay∂x∂ay∂y∂ay∂z∂az∂x∂az∂y∂ax∂z.The matrices G0 and G denote G(t0) and G(t), respectively, and Δt is (t-t0).


## 4. Illustration

STM is not to be used as a propagator and is used between the updates of the state over the duration of orbit determination. However, an experiment of propagation is carried out here over certain duration to ensure that the present model efficiently captures the secular effects. The selected orbit is Molniya orbit witha=26554 kms, e=0.72, and i=63.4 degrees. This orbit is more eccentric and the short periodic effects are predominant. The top dotted line in Figure 2 is the result that depicts the absolute deviation in position (distance) from the STM that includes short and periodic effects given in [7] with respect to (3a) and (3b). The bottom curve is the absolute deviation between the proposed STM (15)–(19) with respect to (3a) and (3b) again as function of time. Both have the same initial conditions. This illustrates that the proposed model is closer to the mean or secular orbital motion given in (3a) and (3b) in the presence of J2 effect neglecting the numerical error due to propagation. The numerical propagation error due to step size is common in both. This validates that the STM derived in Section 3 using the partial derivatives of accelerations in (15) contains secular effects alone.Figure 2
Deviation in position for Molniya orbit.Next a similar exercise as carried out in Figure2 is considering a sun-synchronous polar orbit with a=7168.4, e=0.0011, and i=98.51. The top dashed line in Figure 3 is the result that depicts the deviation in position from the STM that includes short and periodic effects. The bottom curve is the deviation in position from the proposed STM. The illustration confirms that the proposed model is closer to (3a) and (3b).Figure 3
Deviation in position for sun-synchronous orbit.Next example is about illustrating the modeling of (7). The satellite orbit has a=7244.3, e=0.0003, and i=20 degrees. The satellites placed in such orbits are mostly used for metrological purposes. The payload sensor usually covers a wide area as in microwave remote sensing applications. The argument of latitude θ2 corresponding to (3a) and (3b) is first derived. The satellite motion is also obtained as solution of (3a) and (3b) while considering complete perturbations as mentioned in (1). The instantaneous argument of latitude is θa. From (3a) and (3b) the argument of latitude, denoted as θ2, is obtained. The residue (θa-θ2) is plotted against time in Figure 4. It can be seen that this deviation is secular and has large effect along the track.Figure 4
Difference in argument of latitude.The difference between the argument of latitude of the full force modelθa and θ2 over one orbit, that is, 102 minutes in duration, is then fit for the secular effect by a polynomial in a least square sense. For in this case it is(20)θf=θb-0.03624+0.0012t-0.000013t2+0.000000046t3.The same fit is proposed here to extend over successive orbits. However, at the end of each orbital period the ordinate value of difference θb is used for the next orbital period. Here, for example, the value of θb is −0.113 at the start. This approach of orbit model representation can be considered as an alternative to the existing methods for onboard applications [12].The difference between the argument latitude of the full force modelθa and θf in (7) is computed. This residue is illustrated in Figure 5 over four orbits and is observed to be periodic over each orbital period. The initial states for the STM are obtained from the proposed orbit model while using (3a) and (3b) along with the correction in (20). The absolute deviation is with respect to the position of full force secular propagation model. The plot in Figure 6 is similar to that in Figure 1 or Figure 2. The deviation in position in the plot uses the STM in Section 3 and using (15) can be seen to accommodate secular effects in the presence all perturbations.Figure 5
Residue with respect to the polynomial fit.Figure 6
Deviation in position.


## 5. Application

Equations (3a) and (3b) describe the dynamics of the orbiting satellite when oblate earth effect is considered. The relative motion based on geometric approach is given in [9]. They independently make use of the dynamics in (3a) and (3b) for both satellites. The two satellites have identical (a,e,i) in (3a) and distinct (w,Ω,M) in (3b). The relative motion as reproduced from [9] is(21)ΔxΔt=-1+c2im2c2is2cδθ+δΩ+s2im2s2is2cδθ-δΩ+s2im2c2is2c(2θm+δθ+δΩ)+c2im2s2is2c(2θ1+δθ-δΩ)+12simsiscδθ-c2θm+δθ,ΔyΔt=c2im2c2is2s(δθ+δΩ)+s2im2s2is2s(δθ-δΩ)-s2im2c2is2s(2θm+δθ+δΩ)+c2im2s2is2s(2θm+δθ-δΩ)+12simsissδθ+s2θm+δθ,ΔzΔt=-s(im)s(δΩ)c(θs)-sim2ciscδΩ-cimsissθs.Here θm is the instantaneous argument of latitude, that is, θ2 of the main or chief satellite. Similarly θs is the instantaneous argument of latitude of the second or follower satellite and their difference (θm-θs) which is denoted as δθ. The differences in the instantaneous longitude of the ascending node between them is δΩ. The c and s functions denote cos and sin functions and im and is denote the inclinations of the main and follower, respectively. Also Δt denotes the incremental time as used in (3a) and (3b). In (21), Δx,Δy,Δz are relative positions of the second satellite in the orbit frame of the first.Here we shall outline the application of the secular approximationθp from (7) when substituted into (21). This is to enhance (21) to match the secular effects particularly along the track while considering all perturbations. To do so we replace the argument of latitude angle θ2 for both satellites in (21) by an appropriate θp, computed individually using (20) and (7) as in Figure 6, for the main and secondary satellites. These computed θp are then substituted into (21) as θmˇ and θsˇ, respectively, for the main and secondary satellite. Now, (21) is enhanced to match the secular effects particularly along the track of the full force model. It may be noted that in formation flying, (21) is used for guidance and can also be used to derive the initial conditions.


## 6. Conclusion

Approximate orbit model that captures the secular motion in the argument of latitude while considering all perturbation forces is obtained. Approximate state transition matrix is first derived for the orbital motion that matches the secular motion with the effects of oblate earth. This STM can also accommodate the approximate orbit model with all forces. Numerical behavior of the state transition matrix and the orbit models has been provided. Finally a direct application of the orbit model in relative motion of two satellites has been indicated. The suggested approach is simpler for onboard implementation.



---

*Source: 475742-2015-03-10.xml*