Zhang Yonghe,Wang Yamin,Mao Qingyun,Liang Xuwena,

aShanghai Institute of Microsystem and Information Technology,Chinese Academy of Sciences,Shanghai 200050,China

Orbital reference frame estimation with power spectral density constraints for drag-free satellites

Zhang Yonghea,b,c,*,Wang Yaminb,c,Mao Qingyunb,c,Liang Xuwena,b,c

aShanghai Institute of Microsystem and Information Technology,Chinese Academy of Sciences,Shanghai 200050,China

The drag-free satellites are widely used in the field offundamental science as they enable the high-precision measurement in pure gravity fields.This paper investigates the estimation of local orbital reference frame(LORF)for drag-free satellites.An approach,taking account of the combination of the minimum estimation error and power spectral density(PSD)constraint in frequency domain,is proposed.Firstly,the relationship between eigenvalues of estimator and transfer function is built to analyze the suppression and amplification effect on input signals and obtain the eigenvalue range.Secondly,an optimization model for state estimator design with minimum estimation error in time domain and PSD constraint in frequency domain is established.It is solved by the sequential quadratic programming(SQP)algorithm.Finally,the orbital reference frame estimation of low-earth-orbit satellite is taken as an example,and the estimator of minimum variance with PSD constraint is designed and analyzed using the method proposed in this paper.

1.Introduction

Drag-free satellites can eliminate the non-gravitational accelerations,i.e.,solar radiation pressure and atmospheric drag to obtain a free-falling dynamical environment for the floating proof-mass inside the satellite.This kind of satellite is widely used in the measurement of Earth gravity field,the test of Equivalence Principle and the detection of gravitational waves.1Both NASA and ESA have successfully launched drag-free satellites to measure the Earth gravity2–4and proposed the next generation missions,for instance,the NGGM mission5,6and the GRACE follow-on mission.7,8The Microscope mission9for testing the equivalence principle,and LISA mission10,11for detecting gravitational waves have also been scheduled.China has also proposed several space exploration plans based on drag-free technology.12,13

Drag-free concept was originally proposed by Pugh,14and then studied systematically by Lange.15With the development of drag-free missions,a wide variety of studies about the dragfree control has been carried out.16–19The control system con-sists of the inner and outer control loops which take charge of non-gravitational accelerations and attitude pointing,respectively.20The satellite attitude measured by star sensor in earth inertial frame should be aligned with the local orbital reference frame(LORF)which is determined by the satellite position and velocity.The measured attitude is transformed into LORF as the input to the control loop.Due to the stringent requirements on the power spectral density(PSD)of Gravity Field and Steady-State Ocean Circulation Explorer(GOCE)mission,that is,the residual linear non-gravitational accelerations with the PSD below 10-6m/s2/Hz1/2,as well as the angular accelerations with the PSD below 10-6rad/s2/Hz1/2in measurement bandwidth(MBW),the attitude control error of the drag-free satellite is required to reach the level of 10-5rad/Hz1/2in terms of PSD.20Thus,the knowledge of LORF should be estimated at a level better than attitude accuracy in order to minimize the contribution to the overall error budget.21Simultaneously,the high-precision state estimation or accurate measurement of inter-satellites range is the prerequisite of the recovery of gravitational field or the extraction of signals of gravitational waves.Finally,the state estimation is critical for the success of this kind of drag-free satellites and related space science missions.

The satellite position and velocity are measured by onboard GPS,which introduces noises to be filtered.For the state estimation,a series of nonlinear estimators has been proposed in the past.Most of them are nonlinear extensions of the Kalman Filter22which include the Extended Kalman Filter,23Unscented Kalman Filter,24and others.These filters are developed for minimum estimation error and require the knowledge of input noise statistics.However,the drag-free satellite whose main function is the high-precision measurement within certain MBW has special constraints in frequency domain for state estimation.The PSD of estimation error is strictly suppressed below the specification in MBW,but can be increased at the frequency lower than left bound of MBW since the peak attitude is allowed to relax to mrad level.These filters based on optimal estimation error are unable to cover the frequency domain constraints.Thus,a new type of estimator,which has the capability of combining both the PSD constraint in frequency domain and minimum estimation error in time domain,needs to be investigated.

For the attitude estimation,Evers16analyzed the dynamics and estimated the colored measurement noises for GOCE satellite by tuning the covariance matrix of model noises.Canuto25–27proposed the embedded model to deal with the state estimation and drag-free control.For orbital reference estimation,the embedded model consists of high-precision discrete orbital dynamics and disturbed dynamics models.An error feedback gain was designed to construct the time invariant estimator.The PSD constraint of estimation error is met by tuning the eigenvalues of state estimator.Moreover,it doesn’t rely on the statistic knowledge of input error.

In this paper,an extensive study about estimator design is carried out on the basis of the work in Refs.21,24A design method for optimal LORF estimation is proposed to realize the minimum estimation variance with PSD constraints for drag-free satellites.First,the linear time-invariant estimator is established by the time-varying gain matrix.The PSD of input and output signals is connected by the transfer function of the estimator.Second,the PSD constraint is mapped to the requirements on eigenvalues and transfer function by analyzing the magnification and suppression effect on the PSD of input signals.The optimization issue,the minimum estimation error with the PSD constraints,is modeled.It is solved by a local optimal algorithm,that is,the sequential quadratic programming(SQP)algorithm.To validate the adaptability of this method,the LORF of a 250 km low-earth-orbit(LEO)satellite is estimated by the method proposed in this paper.This method is also suitable for the orbital reference estimation of all LEO drag-free satellites.Furthermore,it provides the technology support for the future satellite gravity measurementmission,Space Advanced Gravity Measurement(SAGM),12in China.

2.Background of orbit estimation

2.1.Equations of orbital movement

It is assumed that the drag-free satellite can totally eliminate the non-gravitational accelerations,and then the satellite moves in a pure gravitational field.The acceleration of gravity with J2perturbation is

where J2is set to be 1.08×10-3;R is the Earth equator radius,and R=6378.14 km in this paper;r denotes satellite vector,r= ‖r‖; Γ =diag（1，1，3）; μ is the gravitational constant of the Earth;δg is the residual gravity acceleration;z is the third coordinate component of satellite state.

The orbit of satellite is not circular due to the non-spherical perturbations.Thus,r in Eq.(1)will be time-varying.If the mean orbit height is h,the 1st order expansion of Eq.(1)at the mean radius r=R+h is

where ω0is the mean angular velocity.

Hence,the differential equation of orbit state[r,v]Tcan be written as

where adis the sum of residual non-gravitational accelerations and gravity anomalies.The J2term is included in ∂Ω（r）.

2.2.Discrete-time dynamics

2.2.1.Simplification of discrete dynamics

The time interval of drag-free control is 0.1 s which is far less than the orbital period.Thus,the differential Eq.(3)varies slowly during the 0.1 s.The angular rate within this time interval could be approximated as

where riand ωiare the position vector and angular velocity vector of satellite at time ti,respectively.

For a single coordinate rkand vk(k=1,2,3),the state transition matrix Ackis

where Δt is the time step of discrete dynamic model;Fkis the matrix from Eq.(3),

where ωkiis the angular velocity of the kth coordinate component at time ti.

By replacing the exponential function Eq.(5)with polynomial expansions and taking the first three orders,the discrete dynamic of rkand vkwith units ‘m”could be written as

The above dynamic equations ignore the higher order terms of J2perturbation and the Earth gravity anomalies.A disturbed dynamicmodelisdesigned to compensatethis simplification.

2.2.2.Disturbed dynamics equations

According to the study of Ref.20,there are two kinds of noises in the disturbed dynamics,residual non-gravitational acceleration and Earth gravity anomalies.The PSD analysis of disturbance acceleration of GOCE satellite indicated that the residual non-gravitational acceleration could be approximated as white noise w0k.The main perturbation of Earth gravity is from J2term which has been modeled in the dynamics model in Section 2.2.1.Other non-spherical gravity acceleration could be considered as two random drifts zk=[z1k，z2k]Twhich are accumulative acceleration and are propagated by an integration process.The initial values of random drifts are a pair of white noises w1kand w2kwhich are collected together with w0kinto wk.

Thus,the disturbed dynamics in a single coordinate is

2.3.Closed-loop state estimator

Substitute disturbed Eq.(8)into Eq.(7),and the entire state estimation equation is

and the relevant matrices hold

We denote the vector of GPS measurement with Y（i）and model error e with e=Y（i）-Yk（i）.The relationship between acceleration noises¯wk（i）and the model error e is built by the gain matrix L(i):

In order to estimate the three-dimensional acceleration noises¯wk（i）from the two-dimensional model error e,a 3×2 dimensions gain matrix L(i)is designed.The GPS measurement error,with 30 m(1σ)position error and 0.03 m/s(1σ)velocity error,shows that the position measurement is much noisier than that of velocity.In order to let all measurement noise pass through the feedback channels of noise estimator,the gain matrix L（i）is designed in Ref.21as

Combining Eqs.(9)and(11)yields the feedback state estimator:

2.4.LORF estimation error definition

The illustration of LORF is shown in Fig.1.The inertial frame and LORF are denoted by OIXIYIZIand OOXOYOZO,respectively.

The three axes of LORF are defined by

where r and v are the position and velocity vectors in inertial frame.In Fig.1,~r denotes the vector measured by GPS;^r denotes the estimated vector;and then the LORF estimated error can be expressed by the Euler angles of the transition matrix between frame OOXOYOZOand frame^OO^XO^YO^ZO.ROis defined to denote the real transition matrix from inertial frame to true LORF,and^ROis defined to denote the corresponding estimated transition matrix;and then the relation between ROand^ROcould be expressed as

Fig.1 Illustration of LORF and inertial frame.

3.State estimator with PSD constraint

3.1.Transfer function and PSD

For the linear time-invariant system Eq.(13),the input and output PSDs are linked by the transfer function matrix T（s）,28

where Y（s）is the PSD of input signals,i.e.,measurement values and^Y（s）is the PSD of output signals,i.e.,estimation values.The transfer function matrix of discrete system Eq.(13)is

where s=ei2πfT,with T the time interval of sample,and f the sampling frequency.D is equal to 0 since there is no control force in this dynamic system.For this multi-input and multioutput system Eq.(13),^Y（s）and Y（s）are two-dimensional vector and T（s）is 2× 2 matrix,

Thus, both of position and velocity estimation^Yk（i）=[^rk，^vk]will be influenced by a pair of elements of transfer function T（s）.The relationship is shown as follows:

For the drag-free satellite,such as GOCE,the PSD of orientation error of LORF in the MBW[5×10-3,0.1]Hz is constrained below 4.3 μrad ·Hz-1/2;besides,the variance of entire estimation error is less than 200 μrad in time domain.The illustration of PSD constraint is shown in Fig.2.

As shown in Fig.2,the measurement noise within MBW should be suppressed strictly below 4.3 μrad.However,the PSD constraint in low frequency range is relaxed.Thus,the state estimator in this paper should have two functions:lowpass filter in frequency domain and minimum estimation error in time domain.

For a low-pass filter of linear system,the eigenvalues(pole points)offilter should be within the unit disk of complex plane for stability and are suggested in the neighborhood of(1,0i)for the low-pass character.28The suppression effect is also relevant to both of the zero pole points of transfer function Eq.(18).

3.2.State estimator of minimum error with PSD constraint

For a linear time-invariant system,the state transition matrix A of Eq.(13)should be time-invariant,which is realized by a time-varying gain matrix L（i）with four variables.For a single coordinate,the estimator Eq.(13)is four-dimensional,which is determined by four eigenvalues.That means the gain matrix L（i） can be calculated uniquely by the eigenvalues λi（i=1，2，3，4）of matrix A of the estimator.The calculation procedure is described as follows:(1)calculate the state transition matrix Ackby the orbit estimation values;(2)solve the linear algebraic equation ‖Ack- λI‖ =0 by the given eigenvalues λi（i=1，2，3，4）to acquire the time-varying gain L（i）.

Thus,in respect of time domain,the convergence of closedloop estimator and variance of entire estimation error are dependent on the eigenvalues λ;in respect offrequency domain,the transfer function of system Eq.(13)is determined by the eigenvalues λ and determines the PSD of estimation errors.A design approach for state estimator combining optimal estimation error in time domain and PSD constraint in frequency domain is proposed here.

Fig.2 Illustration of PSD constraint of orbital reference estimation for GOCE satellite.

The approach consists of analysis of PSD constraint and optimization of estimation error.In order to reduce the search space of eigenvalues in optimization process,this approach deals with the PSD constraint by using the transfer function preliminarily.Because the transfer function could evaluate the input PSD analytically without time-domain simulation.The boundaries of eigenvalues are obtained by analyzing the magnification and suppression effects on the input PSD.In particular,the relationship between the starting point of suppression frequency f1and the eigenvalues is obtained.For example,asshown in Fig.2,thePSD constraintis 0.005 Hz≤f≤0.1 Hz and f1should be less than 0.005 Hz.The value off1can be calculated by an ordinary solution algorithm of nonlinear equations.After this calculation procedure,the boundaries for eigenvalues can be determined.From the analysis of eigenvalues and transfer function,the trend of PSD of output signal within a small search space of eigenvalues is obtained.The boundary of eigenvalues,obtained by the trend analysis,could be thought to be located in the neighborhood of optimal eigenvalues and is a good guess for the local optimization algorithm.

Then,the optimization model for minimum estimation error with the eigenvalues constraints is constructed and solved by the local optimization algorithm SQP.During the optimization procedure,the orbit estimation in time domain is executed by Eq.(13).The GPS measurement data are from the highprecision orbital propagation in the dynamic model considering 21×21 WGS84 Earth gravity,Moon,and Sun gravities with Gaussian white noises.The orientation error^eOof LORF is computed by the expression Eq.(16).The PSD of orientation error is computed by the ITPDA Toolbox.29

These four eigenvalues determine both of the zeros and poles of transfer function of system Eq.(18).The successful suppression effect of transfer function on position/velocity probably accompanies the magnification effect on velocity/position in frequency domain,which may lead to the severe PSD of LORF error within MBW and the increase of error variance in time domain.The optimal estimation error in time domain and PSD constraint in frequency domain are two competing factors.Therefore,the eigenvalues should be designed with a compromise between the suppression effect in frequency range and the error variance in time domain.The optimal model can be summarized as follows:

In order to meet the PSD constraint,the optimization model still has the PSD constraint except for the boundary constraints of eigenvalues.The entire design procedure of optimal estimator is shown in Fig.3.

4.Orbital reference estimation for LEO drag-free satellite

4.1.Analysis of eigenvalues margins

Fig.3 Design procedure of optimal estimator with PSD constraint.

The Gaussian white noises are applied to analyze the suppression and magnification effect on input signals.For a 250 km LEO satellite,the noises of GPS for position and velocity are σr=30 m （1σ）and σv=0.03 m/s （1σ）,respectively,and the sample time Ts=1 s,i.e.,sample frequency fs=1 Hz.Then,the measured error of LORF can be evaluated through Eq.(16),and the PSD of error angle is calculated as

Eq.(22)shows that the direct LORF measurement is not compliant with the 4.3 μrad ·Hz-1/2of constraint limitation in MBW,and an estimator is necessary for LORF in real time to improve the reference frame knowledge.

The PSDs of position noise and velocity noise are shown in Fig.4.For a given interval of sampling time,the PSD of white noise is almost distributed equally.The state estimator will suppress the high-frequency signals when the eigenvalues are set in the neighborhood of(1,0i).

Four elements of transfer function T（s）with eigenvalues λi=0.99 （i=1，2，3，4）and 0.999 are separately shown in Fig.5(a)and(b).Its frequency range is[1×10-5,5]Hz.When λi=0.99 （i=1，2，3，4）,the transfer function varying with frequency is shown in Fig.5(a).According to Eq.(20),T1and T2multiply rkand vk,respectively,to generate the PSD of^rk.As seen in Fig.5(a),the amplitudes corresponding to high frequency range,f＞0.01 Hz,are suppressed a lot.However,it must be noted that T4magnifies the amplitude of vkaround the frequency 0.01 Hz.The suppression of LORF error,which is relative to both of the rkand vkerror,is associated with magnification effect in some frequency range.The starting point of suppression frequency f1is nearly equal to 0.01 Hz.In order to expand the suppression frequency ranges to lower frequency,the eigenvalues should be closer to(1,0i).

Fig.4 PSDs of Gaussian white noises for position and velocity noises.

Fig.5(b)shows the transfer function with λi=0.999 （i=1，2，3，4）.The starting point of suppression frequency f1shifts to about 0.001 Hz.It satisfies the requirement of GOCE MBW.However,the T4of transfer function magnifies the amplitude of vkaround the frequency 0.001 Hz much more.Compared to Fig.5(a),the magnification ratio of vkis greater,which will generate worse state estimation in time domain.

The PSDs of estimation error for LORF with the eigenvalues of 0.99 and 0.999 are shown in Fig.6(a)and(b).The frequency range of peak value and suppression frequency range agree with the above analysis of transfer function.When λi=0.99 （i=1，2，3，4）,the PSDs of y and z error of LORF are suppressed from 0.01 Hz to high frequency and rise at the neighborhood of 0.01 Hz.When the eigenvalues increase to 0.999,the suppression frequency range is expanded to the low frequency 0.001 Hz.However,the amplitudes of low frequency range around 0.001 Hz will rise up,which results from the magnification effect of transfer function on the velocity error.

From the analysisabove,we can find that greater eigenvalues produce wider suppression frequency ranges.However,the increase in eigenvalue magnitude will magnify the PSD of input signal in the middle frequency range.It may break the mission constraints in frequency and time domains.For a LEO satellite,such as GOCE and SGAM,greater than 0.99 and less than 1 would be a good choice for the eigenvalues.

Fig.5 Magnitude of transfer function in frequency range[2×10-5,1]Hz.

4.2.Optimization design of state estimator

This section will optimize the estimator from the view of minimum variance of estimation error and the view of minimum variance with PSD constraint by using the approach proposed in Section 3.The orbit of GOCE is taken as an example.Its orbit is a near circular one with orbital altitude 250 km.In the numerical simulation,the initial eccentricity is 0,which is not constant and will evolve periodically in the dynamical model with J2perturbation.The constraints are as follows:(1)MBW 0.005 Hz≤f≤0.1 Hz;(2)amplitudes of PSD within the MBW less than 4.3 μrad ·Hz-1/2;(3)the error variance in the entire time domain less than 200 μrad.

4.2.1.Minimum variance of estimation error

By using the approach proposed in Section 3 without the PSD constraint,the minimum variance of estimation error in the time domain is obtained.Twelve eigenvalues for the rkand vkof three coordinates are{0.995483,0.996608,0.991828,0.999871,0.999736,0.992393,0.998912,0.994079,0.997859,0.999900,0.990314,0.990594}.The errors of r1and v1of inertial frame are shown in Fig.7.It is clear that compared to the GPS errors,the state estimation errors are almost decreased by one order of magnitude.

Fig.6 PSDs of LORF estimation error with eigenvalues of 0.99 and 0.999.

By converting the state errors of inertial frame to the LORF,the errors of LORF are obtained.The optimal error variances for three coordinates are less than 0.8 μrad.Their PSDs are shown in Fig.8.The MBW is marked.Some of the error amplitudes associated with the frequencies between 0.005 Hz and 0.01 Hz are greater than the required upper limit 4.3 μrad.It means that greater eigenvalues are required to expand the suppression frequency range to a lower frequency.

Thus,the estimator of minimum error variance may not meet the PSD constraint in the frequency domain.

4.2.2.Minimum variance with PSD constraint

By using the optimization method proposed in Section 3,the minimum variance of estimation error with PSD constraint is designed. Twelve eigenvalues are {0.998981, 0.999923,0.996176,0.997671,0.996533,0.999264,0.997781,0.997053,0.997544,0.999030,0.996703,0.999105}.Obviously,the magnitude of eigenvalues increases slightly in comparison with that of minimum variance estimator.The transfer function corresponding to the optimal eigenvalues is shown in Fig.9.We can see the starting point of suppression frequency is about 0.0002 Hz.The suppression frequency range includes the MBWofGOCE,0.005 Hz≤f≤0.1 Hz.Themaximum amplitude associated with the MBW is 3.7 μrad which meets the PSD constraint.Compared to the transfer function associated with λi=0.999 （i=1，2，3，4）which is shown in Fig.5(b),the magnification ratio of velocity estimation is reduced a lot.The optimal state estimator can reduce the amplitude of PSD of estimation errors within the MBW and meet the PSD constraint.The PSDs of LORF error are shown in Fig.10.

Fig.7 GPS noises and estimation errors of position and velocity in inertial frame.

Fig.8 PSDs of minimum estimation error of LORF.

Fig.9 Transfer function of optimal estimator with PSD constraint.

Fig.10 PSDs of minimum LORF errors with PSD constraint.

Due to the high-pass characteristic of transfer function,almost all of the low frequency signals of GPS noises pass the estimator.Thus,the amplitude associated with low frequency range in Fig.10 is big.This phenomenon is also validated by the errors of time-domain x and y components of inertial frame shown in Fig.11.We can find that the lowfrequency noises of GPS remain after filtered by the optimal estimator.

5.Conclusions

The orbital reference frame estimation issue for the drag-free satellite is investigated in this paper.A design method for optimization estimator is proposed to deal with the orbital reference frame estimation with PSD constraint in frequency domain.Some conclusions drawn from the results are summarized as follows:

(1)Bigger eigenvalues will result in better suppression effect on position estimation but higher magnification effect on velocity estimation.From the analysis of transfer function,the eigenvalue ranges are suggested to be[0.99,1].

(2)Estimator with minimum estimation error may not meet the PSD constraint.The procedure for estimator design proposed in this paper can optimize the eigenvalues to reconcile the suppression and magnification effects and meet the PSD constraint.

Fig.11 Optimal estimation errors of position and velocity in inertial frame.

(3)The method proposed in this paper could also be applied to other next generation LEO drag-free satellites in the future.

Acknowledgments

This study was co-supported by the Open Fund of Joint Key Laboratory of Microsatellite of CAS(No.KFKT15SYS1)and the Innovation Foundation of CAS(No.CXJJ-14-Q52).

1.DeBra DB.Drag-free control for fundamental physics missions.Adv Space Res 2003;32(7):1221–6.

2.Nkwater MR,Floberghagen R,Haagmans R,Muzi D,Popescu A.GOCE:ESA’s first Earth explorer core mission.Space Sci Rev 2003;108(1–2):419–32.

3.Tapley BD,Bettadpur S,Ries JC,Thompson PF,Watkins MM.GRACE measurements of mass variability in the Earth system.Science 2004;305(5683):503–5.

4.Reigber C,Lu¨hr H,Schwintzer P.CHAMP mission status.Adv Space Res 2002;30(2):129–34.

5.Cesare S,Mottini S,Musso F,Parisch M,Sechi G,Canuto E,et al.Satellite formation for a next generation gravimetry mission.In:Sandau R,Roeser H,Valenzuela A,editors.Small satellite missions for earth observation.Berlin:Springer;2010.p.125–33.

6.Panet I,Flury J,Biancale R,Gruber T,Johannessen J,van den Broeke MR,et al.Earth system mass transport mission(e.motion):A concept for future earth gravity field measurements from space.Surv Geophys 2013;34(2):141–63.

7.Zheng W,Xu HZ,Zhong M,Yun MJ.Requirements analysis for future satellite gravity mission improved-GRACE.Surv Geophys 2015;36(1):87–109.

8.Loomis BD,Nerem RS,Luthcke SB.Simulation study of a follow-ongravitymissiontoGRACE.JGeodesy2012;86(5):319–35.

9.Touboul P,Me´tris G,Lebat V,Robert A.The MICROSCOPE experiment,ready for the in-orbit test of the equivalence principle.Class Quant Grav 2012;29(18):1–6.

10.Danzmann K.LISA–An ESA cornerstone mission for the detection and observation of gravitational waves.Adv Space Res 2003;32(7):1233–42.

11.Fichter W,Gath P,Vitale S.LISA pathfinder drag-free control and system implications.Class Quant Grav 2005;22:139–48.

12.Hu WR,Xu HZ.GRACE fellow-on in China—space advanced gravity measurement program.Saarbru¨cken:Lambert Academic Publishing;2015.p.23–38.

13.Luo L,Chen LS,Duan HZ,Gong YG,Ji JH,Liu Q,et al.TianQin:A space-borne gravitational wave detector.Class Quant Grav 2016;33(3):035010.

14.Pugh G.Proposal for a satellite test of the coriolis prediction of general relativity.Washington,D.C.:Weapons Systems Evaluation Group Research Memorandum,NASA;1959.Report No:NASA-No.11.

15.Lange B.The control and use of drag-free satellites[dissertation].California:Stanford University;1964.

16.Evers W.GOCE dynamical analysis and drag free mode control.ESA/ESTEC;2001.Report No:DCT-2004.44.

17.Theil S,From M,Pettazzi L,Scheithauer S.Formations of dragfree satellites.Proceedings of the 56th international astronautical congress;2005 Oct 17-21;Fakuoka,Japan;2005.

18.Zhang YH,Liang XW,Zhang J,Hu QL.Finite time relative position control for drag-free dual-satellite serial formation.J Astronaut 2015;36(8):923–31[Chinese].

19.Wu SF,Fertin D.Spacecraft drag-free attitude control system design with quantitative feedback theory.Acta Astronaut 2008;62(12):668–82.

20.Canuto E,Massotti L.All-propulsion design of the drag-free and attitude control of the European satellite GOCE.Acta Astronaut 2009;64(2–3):325–44.

21.Canuto E,Massotti L.Local orbital frame predictor for LEO drag-free satellite.Acta Astronaut 2010;66(3–4):446–54.

22.Kalman RE.A new approach to linear filtering and prediction problems.J Basic Eng 1960;82:167–79.

23.Gelb A.Applied optimal estimation.Cambridge:MIT Press;1974.p.182–203.

24.Wan EA,Rudolph VDM.The unscented Kalman filter for nonlinear estimation.IEEE proceedings of the symposium 2000 on adaptive systems for signal processing,communication and control;2000 Oct 1-4;Lake Louise,Alberta.Piscataway(NJ):IEEE Press;2000.

25.Canuto E.Embedded model control:outline of the theory.ISA Trans 2007;46(3):363–77.

26.Canuto E,Martella P,Sechi G.Attitude and drag control:an application to the GOCE satellite.Space Sci Rev 2003;108(1–2):357–66.

27.Canuto E.Drag-free and attitude control for the GOCE satellite.Automatica 2008;44(7):1766–80.

28.Proakis JG,Manolakis DG.Digital signal processing–principles,algorithms,and applications.4th ed.New York:Pearson Education;2007.p.243–7.

29.LTPDA:A MATLAB toolbox for accountable and reproducible data analysis[Internet].Potsdam:Max Planck Institute for Gravitational Physics;2015[cited 2016 Apr 5].Available from:www.elisascience.org/ltpda.

Zhang Yonghe received the B.S.and M.S.degrees from Harbin Institute of Technology in 2001 and 2003,respectively.He is currently an associate professor at Shanghai Engineering Center for Microsatellites.His current research interests are drag-free control and satellite system design.

bShanghai Engineering Centre for Microsatellites,Shanghai 201210,China

cJoint Key Laboratory of Microsatellite,Chinese Academy of Sciences,Shanghai 201210,China

8 January 2016;revised 16 February 2016;accepted 24 March 2016

Available online 21 October 2016

Drag-free satellite;

Orbital reference frame;

Power spectral density;

State estimation;

Transfer function

Ⓒ2016 Production and hosting by Elsevier Ltd.on behalf of Chinese Society of Aeronautics and Astronautics.This is an open access article under the CC BY-NC-ND license(http://creativecommons.org/licenses/by-nc-nd/4.0/).

*Corresponding author.Tel.:+86 21 50735181.

E-mail address:yonghe.zhang@mail.sim.ac.cn(Y.Zhang).

Peer review under responsibility of Editorial Committee of CJA.

Production and hosting by Elsevier

http://dx.doi.org/10.1016/j.cja.2016.10.013

1000-9361Ⓒ2016 Production and hosting by Elsevier Ltd.on behalf of Chinese Society of Aeronautics and Astronautics.

This is an open access article under the CC BY-NC-ND license(http://creativecommons.org/licenses/by-nc-nd/4.0/).