Modeling and Prediction of Momentum Wheel Speed Data

To solve the problems of data loss and unequal interval of momentum wheel (MW) speed during a satellite stable operation, this paper presents a multidimensional AR model. A Lagrange interpolation method is used to convert measurements to equal interval data, and the FFT algorithm is adopted to calculate the period of MW speed variation. The long data sequence is converted into multidimensional time series, based on the equal interval data and the period. A multidimensional AR model is established, and the least square method is used to estimate the model parameters. The future data trend is predicted by the proposed model. Simulation results show that the prediction algorithm can achieve the across cycle prediction of the MW speed data.


Introduction
The momentum wheel (MW) is the key actuator of the Satellite Attitude Control System (SACS), which includes bearing assembly, motor, casing, and wheels. In the long running of the three-axis stabilized GEO satellite, the MW speed data has a periodical change characteristic. Through long-term monitoring for MW speed of a three-axis stabilized GEO satellite, [1] found that its variation can be sorted into the shortperiod terms (in days) and the long-period terms (in years), and the periodic variation of the MV is also analyzed. The equation of the MW speed is obtained by analyzing the operation principle of the bias MW, and the graphic of MW speed change is described in [2]. In [3], the authors considered that the reliability information of MW can be provided by collecting degradation data when there exists certain performance characteristics that degrade over time; they develop a reliability modeling and life estimation approach for MW used in satellites based on the expectation maximization (EM) algorithm from a Wiener degradation model. GmbH has installed 588 MWs in 235 satellites, which is required fault-free operation more than 2200 years [4].
Under the complex working environment, such as water mist, salt spray, corrosion, and aging resulting from fuels containing sodium, sulfur, vanadium, etc. which are burned, the dynamic characteristics of speed sensors are severely influenced, which causes the problems that the data detected is inaccurate or could not be detected in a certain time that means false detection or data loss [5]. Furthermore, the data detected may be equal interval incompletely because of the unstable precision of a sensor. As we know, the rotational speed is an important state parameter of the momentum wheel. Once the speed sensors occur to faults, it will cause a wrong operation of control system and even lead to accidents. The Telstar satellite of America was lost early on 21 February 1963 because of the abnormal data of a sensor caused by cosmic rays. This was the first spacecraft loss due to radiation effects [6].
The remainder of this paper is organized as follows. In Section 2, we briefly analyze the change characteristic of the MW speed data, using Lagrange interpolation and Fast Fourier Transform Algorithm (FFT) to determine the period of MW speed data. Section 3 establishes the multidimensional model to predict the MW speed data and verifies the validity of the model. Finally, the paper concludes with some brief, summary remarks for the present work.

Obtainment of the Variation Period of the
Momentum Wheel Date communication satellite is time series data. Assuming that the observed MW speed is yðt i Þ at time t i , written as y i , the change curve of MW speed is shown in Figure 1. Because of accidental factors such as cosmic ray and accidental contact failure of circuit, the sensor data is abnormal [7]; for example, the observed data is incompletely equal interval distribution and the partially observed data loss. Meanwhile, it can be seen that the MW speed is approximately periodic and its amplitude grows linearly. In order to predict MW speed, the Lagrange interpolation method [8] is used to calculate the MW speed, so that the unequal interval data is converted into equal interval data.

Preprocessing of Equal Interval Interpolation of Momentum
Wheel Speed Data. The time series data of the observed MW speed is Y = ½ðt 1 , y 1 Þ, ðt 2 , y 2 Þ, ⋯, ðt N , y N Þ T , where N is the total number of data. However, Y is the unequal interval time series data, and the forecast model based on equal interval time series data. So Y has to convert into equal interval time series data by the interpolation method. Firstly, set the interpolation interval p. To interpolate the data segment from ðt 1 , y 1 Þ to ðt m+1 , y m+1 Þ in Y, the Lagrange interpolation formula [9] can be given by where t is the interpolation time; t 1 , t 2 , ⋯, t m+1 are the observation time which satisfy t 1 < t 2 <⋯<t m+1 . According to formula (1), it yields where t 1 ′ and t m+1 ′ are the beginning time and the ending time after interpolation, respectively; the value of int ððt m+1 − t 1 Þ * p −1 Þ is an integer multiple relationship between the difference of time ðt m+1 − t 1 Þ and the interpolation interval p. So the time range after interpolation is Similarly, set the same interpolation interval to interpolate the second segment data from ðt m+1 , y m+1 Þ to ðt 2m+1 , y 2m+1 Þ in Y. The time range of the original data is ðt m+1 , t 2m+1 Þ.
Step 2. The interpolation time data based on formulas (2) and (3), and it is equal.
Step 3. Lagrange Interpolation method is used to interpolate the MW speed data by formula (1).
Step 4. When i > 1, the observed data is y Step 5. The interpolation time data based on formulas (6) and (7), and it is also equal.
Step 6. The Lagrange interpolation method is also used to interpolate the MW speed data by formula (1).
The flow chart of the equal interval interpolation method is shown in Figure 2. 2.3. Determination of the Period of the Momentum Wheel Speed Data. The equal interval data can be calculated by the Lagrange interpolation method, and it is approximately periodic. In order to determine the period of the MW speed data, the Fast Fourier Transformation (FFT) method is adopted to calculate the frequency of the equal interval data Y ′ = ½ðt 1 ′, y 1 ′Þ, ðt 2 ′, y 2 ′Þ, ⋯, ðt N ′ , y N ′ Þ T . If the calculated frequency is a constant, it shows that the MW speed data has the characteristic of approximate periodicity. Setting the time interva1 p = 0:1 s, then the sampling frequency is f s = 1/p = 10 Hz. The FFT method [10] is used to analyze the frequency of the MW speed data. The frequency diagram of MW data is shown in Figure 3.
According to the analysis of data and Figure 3, the amplitude of the 35th data is maximum, so the frequency of MW speed data f = 0:0416 Hz and the period of the MW speed T = f −1 = 24:02 s ≈ 24 s.

Modeling and Prediction of Momentum
Wheel Data E AIC = ln δ 2 + 2 i + j n , where σ 2 is the variance of error between the model and the observed data, i and j are the order of the AR and MA model, respectively, and n is the length of data sequence for modeling and satisfies n = 5 * i + j + 1.
The loss value of AIC and BIC criterion under the length n of the data sequence is shown in Figure 4. We can see that the loss function values of AIC and BIC are relatively smoother between n = 10 and n = 18, which means the robustness of the ARMA model is better. When n = 16, the loss values of AIC and BIC are minimum. Combining n = 5 * i + j + 1 with i > j, it can be calculated when i = 3 and j = 0, and the following study will be based on the third order AR model.

AR Modeling of the Momentum Wheel Speed Data.
In order to extract the MW speed data conveniently, we convert the equal interval data sequence Y ′ = ½ðt 1 ′ , y 1 ′ Þ, ðt 2 ′ , y 2 ′ Þ, ⋯, ðt N ′ , y N ′ Þ T into Y ″ = ½ð1, y 1 ′ Þ, ð2, y 2 ′ Þ, ⋯, ðN, y N ′ Þ T , and the period of Y ″ is calculated as T ′ = T/p = 240. The curve of MW speed data is similar with the sinusoidal function. In order to reduce the calculation, the data of 0, quarter, a half, and three quarter periods are selected as characteristic points expressed by h 1 , h 60 , h 120 , and h 180 . We find that each of them is basically linear in Figure 5. Therefore, an AR model can be established for each element [12][13][14]. As discussed above, the loss function is minimized with a third order AR model, the AR equation is shown as follows: where a 1 , a 2 , and a 3 are the AR model coefficients, feðkÞg is a white Gaussian noise with zero mean and variance δ e 2 = 0:09, Assume φ T ðkÞ = ½xðk − 1Þ, xðk − 2Þ, xðk − 3Þ and θ = ½a 1 , a 2 , a 3 T , the formula (9) can be written as a vector form by According to (10), we can obtain Data segments that need to be interpolated (t j ,y j ), (j = i,i+1,...,i+m) Data segments that need to be interpolated (t j ,y j ), (j = i,i+1,...,i+m)

International Journal of Aerospace Engineering
The least-square method is used to calculate the minimum of RSS (residual sum of squares) for unknown parameter vector θ. The performance function JðθÞ is The purpose of equation (12) is to obtain the estimation value b θ, when the function JðθÞ is minimized.
Formula (10) is rewritten as a vector form , The least-square estimation of θ [15] at time k is 3.3. Verification of the MW Speed AR Model. The first 6000 data is selected as a training dataset, and others is as a test data; the estimated value of model coefficient can be obtained by equation (15).
Based on formulas (10) and (15), the value of xðkÞ is obtained which can be used to predict the value of xðk + 1Þ. Therefore, the MW speed data in the next period can be obtained by the AR model. The forecast model is shown as follows:x The predicted value and the measured value are shown in Figure 6. It can be seen that the predicted data is in accordance with the measured data. The error between the measured value and the predicted value is shown in Figure 7. The results show the effectiveness and feasibility of the proposed algorithm.
The standard deviation of the residuals is used to verify the accuracy of the forecast model, and it is expressed by where x i is the actual observed data, andx i is the predicted data of the forecast model. The results are shown in Table 1.