Time-varying state-space model identification of an on-orbit rigid-flexible coupling spacecraft using an improved predictor-based recursive subspace algorithm

Acta Astronautica 163 (2019) 157–167

Contents lists available at ScienceDirect

Acta Astronautica journal homepage: www.elsevier.com/locate/actaastro

Time-varying state-space model identification of an on-orbit rigid-flexible coupling spacecraft using an improved predictor-based recursive subspace algorithm

T

Zhiyu Nia, Jinguo Liub,∗, Shunan Wuc, Zhigang Wuc a

College of Aerospace Engineering, Shenyang Aerospace University, Shenyang, China State Key Laboratory of Robotics, Shenyang Institute of Automation, Chinese Academy of Sciences, Shenyang, China c School of Aeronautics and Astronautics, Dalian University of Technology, Dalian, China b

A R T I C LE I N FO

A B S T R A C T

Keywords: Parameter identification State-space model Recursive subspace method Time-varying system Affine projection sign algorithm State estimation

Spacecraft control problems frequently require the latest model parameters to provide timely updates to the controller parameters. This study investigates the recursive identification problem in a the time-varying statespace model of an on-orbit rigid-flexible coupling spacecraft. An improved recursive predictor-based subspace identification (RPBSID) method is presented to increase on-orbit identification efficiency. Compared with the classical RPBSID and other subspace methods, the improved RPBSID applies the affine projection sign algorithm. Accordingly, the system state variables can be determined directly via recursive computation. Thus, the proposed algorithm does not require constructing the corresponding Hankel matrix or implementing singular value decomposition (SVD) at each time instant. Consequently, the amount of data used in the identification process is reduced, and the computational complexity of the original method is decreased. The time-varying state-space model of the spacecraft is estimated through numerical simulations using the classical RPBSID, improved RPBSID, and SVD-based approaches. The computational efficiency and accuracy of the three methods are compared for different system orders. Computed results of the test response demonstrate that the improved RPBSID algorithm not only achieves sufficient identification accuracy but also exhibits better computational efficiency than the classical methods in identifying the parameters of the spacecraft time-varying state-space model.

1. Introduction The model parameters of spacecraft have been identified in several studies [1–6], including the on-orbit experiments on the system statespace model of the Galileo spacecraft and the Hubble space telescope [7,8]. The identified model parameters not only enable tracking and monitoring of the on-orbit working conditions of spacecraft but can also be used to provide a reference for the design and correction of system controllers [9,10]. A spacecraft typically has large and flexible appendages, and thus, the corresponding system modeling, identification, and control problems have become increasingly complex [11–15]. However, on-orbit experiments have been mostly implemented based on a linear time-invariant (LTI) system. In actual operations, the system state-space model parameters of spacecraft may vary with time due to structural changes during on-orbit operation, such as the rotation and deployment of appendages [16,17], docking with a satellite [18], or capturing another moving body. Therefore, we aim to implement an on∗

orbit identification experiment without affecting the normal working conditions of spacecraft [19]. Accordingly, the model parameters of onorbit spacecraft must be accurately identified for linear time-variant (LTV) systems. The eigensystem realization algorithm (ERA) and subspace model identification (SMI) method have been successfully applied to estimate the LTI state-space model of various spacecraft in available on-orbit identification methods for spacecraft model parameters [3–8,20–22]. These approaches construct a Hankel matrix using complete input and output (IeO) data sequences over a certain period. In addition, singular value decomposition (SVD) is implemented at each time step to extract the observability matrix. The system state-space model matrices can then be obtained. The classical ERA and SMI are primarily based on LTI systems; for LTV systems, repeated experimental methods, such as the pseudo-modal subspace method and the time-varying ERA (TV-ERA) [23,24], are frequently used to identify the time-varying system model [25,26]. In repeated experimental methods, multiple experiments are

Corresponding author. E-mail address: [email protected] (J. Liu).

https://doi.org/10.1016/j.actaastro.2018.11.008 Received 31 July 2018; Received in revised form 30 September 2018; Accepted 8 November 2018 Available online 09 November 2018 0094-5765/ © 2018 IAA. Published by Elsevier Ltd. All rights reserved.

Acta Astronautica 163 (2019) 157–167

Z. Ni et al.

algorithm. In addition, the computational efficiency of the proposed method in identification for different system orders is discussed. In Section 5, several conclusions drawn from the study are finally presented.

performed during identification and the system undergoes the same time-varying change for different inputs during each experiment. Then, the LTV state-space model parameters are determined. However, these identification methods are based on SVD (hereafter referred to as SVDbased methods), and their procedures require considerable computational costs. If the system order is high, then the computational complexity for SVD is huge. Certain control problems, such as self-adaptive control, often require obtaining the latest model parameters to update the controller parameters in real time. Frequently used identification algorithms that are based on SVD are unsuitable for the online updating of controller parameters. Consequently, a faster model identification method for on-orbit spacecraft must be developed. A series of recursive identification algorithms based on theory of signal subspace projection has been proposed and used to identify the time-varying model parameters of spacecraft [27–30]. A new model estimation method, called recursive predictor-based subspace identification (RPBSID), was developed by Houtzager et al. to improve the noise immunity of existing recursive subspace methods [31]. In this recursive algorithm, a vector autoregressive predictor is adopted to provide an asymptotically consistent estimate. Furthermore, adaptive filtering and recursive least squares are used to obtain the time-varying state-space model parameters in real time and to overcome the influence of noise. Therefore, resistance to noise is increased. Although the RPBSID method has been already used in recent years to identify system model parameters in certain aeronautical flexible-wing structures [32–34], few studies have focused on identifying spacecraft state-space model parameters. In addition, the RPBSID method should construct a generalized Hankel matrix at each time instant to compute the system state vector; thus, a huge amount of data is required. To improve the computational efficiency of the original algorithm and maintain satisfactory identification accuracy, an improved RPBSID method is developed in the current study, and the corresponding recursive identification of the parameters of the spacecraft time-varying state-space model is investigated. A rigid-flexible coupling dynamics model of the system is then established considering the coupling influence between spacecraft attitude motion and flexible structure vibration. To improve computational efficiency, recursive least squares technology and the affine projection sign algorithm (APSA) are applied to solve the system state vector recursively without constructing the corresponding generalized Hankel matrix at each time instant. Finally, the system time-varying state-space model parameters are identified using the system's IeO signals. The amount of data required in the identification process is reduced, and the computational complexity of the original method is decreased using the improved RPBSID method. In the simulations, two on-orbit cases that may cause the model parameters to vary with time are considered, and the time-varying statespace model parameters of a satellite are estimated by using the classical RPBSID, improved RPBSID, and SVD-based approaches. Furthermore, the same test inputs are applied to the original and identified state-space models under the initial state conditions. The system responses with respect to the test inputs of the three methods are then compared to verify the accuracy of the identified model parameters. The computational results demonstrate that the improved RPBSID algorithm not only achieve sufficient identification accuracy, but also higher computational efficiency than those of the classical RPBSID and SVD-based series methods. The proposed method can be used to identify the spacecraft time-varying state-space model parameters. The remainder of this paper is organized as follows. In Section 2, the dynamical modeling of a rigid-flexible coupling spacecraft is reviewed, and the time-varying state-space equation of the spacecraft is provided. In Section 3, an improved RPBSID algorithm is proposed, and the procedure for identifying the time-varying model parameters using the method is summarized. In Section 4, the system's IeO signals that are required for the identification are designed, and two on-orbit cases are simulated to validate the identification ability of the improved RPBSID

2. Modeling description of rigid-flexible coupling spacecraft The structures of rigid-flexible coupling spacecraft can be constructed in the form of a central rigid body with N large flexible appendages, such as solar panels or antennas. We define ψ = [ψx , ψy , ψz ]T ∈ 3 × 1(please move the math to the next paragraph) and ω = [ωx , ω y , ωz ]T ∈ 3 × 1 as the attitude angle and angular velocity vectors of the spacecraft, respectively, where subscripts x, y, and z denote the variables for the roll, pitch, and yaw directions, respectively. Then, the following assumptions are made in the study. (1) The origin of the body frame is located at the mass center of the entire spacecraft, thereby implying that the spacecraft translation motion can be decoupled from its rotation motion. (2) The attitude rotation motion of the spacecraft is slow, and its angular velocity is small, i.e., ω ≈ ψ˙ . On the basis of the slow-motion assumption, several nonlinear terms in equation derivation will be considered as high-order small quantities, and thus, are disregarded. (3) This study does not consider the effects of orbital factors and gravity gradients on spacecraft structural identification. Therefore, the kinetic energy Ttotal and potential energy Vtotal of the system can be written as

Ttotal =

Vtotal =

1 T 1 T s˙ ms˙ + ψ˙ Jψ˙ + 2 2

1 2

N

N

∑ s˙ TPi η˙ i + ∑ ψ˙TFi η˙i + i=1

i=1

1 2

N

∑ η˙ iT η˙i i=1

(1)

N

∑ η˙iT Ωi2 η˙˙i

(2)

i=1

3 × 3

is the inertia where m is the total spacecraft mass, and J (t ) ∈ matrix. s = [sx , s y, sz ]T ∈ 3 × 1 is the translational displacement vector of the spacecraft. ηi ∈ κi× 1 and Ωi2 ∈ κi× κi (i = {1,2, ..., N } ) are the modal coordinate and modal stiffness matrix of the ith appendage, respectively, where parameter κi denotes the selected number of vibrational modes of the ith appendage. P ∈ 3 × κi and F ∈ 3 × κi denote the modal momentum and modal angular momentum influence coefficient matrices, respectively. Thus, the Lagrange function of this system can be expressed as

LF = Ttotal − Vtotal

(3)

If a generalized coordinate vector X¯ is defined as T X¯ = [ s T ψT η1T η2T ⋯ ηNT ]

(4)

then Eq. (3) is substituted into the following Lagrange equation:

∂LF d ⎛ ∂LF ⎞ = Q¯ ⎜ ⎟ − ∂X¯ dt ⎝ ∂X¯˙ ⎠

(5)

where Q¯ is the system generalized force, and the dynamics equations of the rigid-flexible coupling spacecraft can be obtained as [35,36] N

ms¨ +

∑ Pi η¨i = f (t ) i=1

(6)

N

Jψ¨ +

∑ Fi η¨i = u (t ) i=1

FiT ψ¨ + η¨i + 2ζi Ωi η˙ i + Ωi2 ηi = 0, i = {1,2, ..., N }

(7) (8)

where f (t ) and u (t ) are the control force and torque, respectively; and ζi is the corresponding damping ratio of the ith appendage. Notably, this study does not consider the translational motion of the spacecraft but is concerned only with its rotational motion. Moreover, the angular velocity of the spacecraft is assumed to be extremely small. If the 158

Acta Astronautica 163 (2019) 157–167

Z. Ni et al.

xk + 1 = Ak xk + Bk uk + Kk ek

spacecraft structural configuration changes during on-orbit operation, then corresponding variations in the system dynamics model parameters will occur, and matrices J and F in Eqs. (7) and (8) will vary with time. Thus, the rigid-flexible coupling equation of the spacecraft with N appendages is given by

J (t ) ψ¨ +

∑ Fi (t ) η¨i = u (t )

yk = Ck xk + ek

n × m

¨ + η¨ + 2ζ Ωi η˙ + Ωi2 η = 0, i = {1,2, ⋯, N } FiT (t ) ψ i i i i

(10)

The following vector δ is selected as T

δ = [ ψT η1T η2T ⋯ ηNT ]

~ ~ xk + 1 = Ak xk + Bk z k

(18)

yk = Ck xk + ek

(19)

~ ~ where Ak = Ak − Kk Ck , Bk = [Bk , Kk ], and z k = [ukT, ykT ]T . A VARX predictor is defined as follows:

(11)

Eqs. (9) and (10) are rewritten as

p

M (t ) δ¨ + Eδ˙ + Kδ = Lu (t )

yk =

(12)

J (t )

F1 (t ) ⋯ Iκ1× κ1 ⋱

Fi (t )

⎤ ⎥, ⎥ IκN × κN ⎦

⎡ 03 × 3 ⎤ ⎢ ⎥ Ω12 =⎢ ⎥, L = ⋱ ⎢ ⎥ 2 ΩN ⎥ ⎢ ⎣ ⎦

where the Markov matrices described by

⎡ 03 × 3 ⎤ ⎢ ⎥ 2ζ1 Ω1 E=⎢ ⎥, K ⋱ ⎢ ⎥ 2ζN ΩN ⎦ ⎣

Ξk(M) −i

I3 × 3 ⎡ 0κ × 3 ⎤ ⎢ ⋮1 ⎥ ⎢ 0κ × 3 ⎥ N ⎦

(13)

0 0 I ⎤ ⎡ ⎤ Ac (t ) = ⎡ −1 (t ) K − M−1 (t ) E , Bc (t ) = M−1 (t ) L − M ⎣ ⎦ ⎣ ⎦

if i = 2,3,⋯, p (21)

(22)

Ξk = Ξk − 1 + (yk − Ξk − 1 φk )(φkT Zk )

(23)

Zk = (Zk − 1 − Zk − 1 φk (β1 I + φkT Zk − 1 φk )−1φkT Zk − 1)/ β1

(24)

where β1 is a forgetting factor that satisfies 0 ≪ β1 ≤ 1. The initial value of the matrix Zk is selected as Z0 = (1/ α1 ) Ip (m + r ) , where α1 > 0 is the initial regularization value to ensure that the problem is well conditioned in the first iterations. The value of the VARX parameter matrix Ξk can be determined using Eq. (23). In the next section, the system state vector xk is estimated recursively using the VARX parameter Ξk .

I is a unit matrix. If the attitude signal of the spacecraft and the vibration displacement signal of the appendages are selected as the output signal y (t ) ∈ m × 1, where m is the number of system output, then the measurement equation of the spacecraft is written as (14)

3.2. Estimation of system state vector using APSA

is the output matrix. The matrix H is

I ⎤, i = {1,2, ⋯, N } H=⎡ ⎣ Φi ⎦

if i = 1

T

where x (t ) = [ δ T (t ) δ˙ T (t ) ] ∈ n × 1 is a state vector, and parameter n is the system order of the state-space model. Hence, the system matrix Ac (t ) ∈ n × n and input matrix Bc (t ) ∈ n × r are expressed as follows:

where C ∈

( i = {1,2,⋯,p − 1, p} ) are

T T T (M) (M) where Ξk = [Ξk(M) − p , Ξk − p + 1, ⋯, Ξk − 1 ] and φk = [ z k − p z k − p + 1 ⋯ z k − 1 ] . There, by using the adaptive filtering technique [37], the least squares recursive form of matrix Ξk ∈ m × p (m + r ) at each time instant is obtained as

T

m × n

~ Ck Bk − 1 ~ ~ ~ = C A ⋯A k − i + 1 Bk − i ⎨ k k−1 ⏟ i−1 ⎩ ⎧

∈

m × (m + r )

yk = Ξk φk

Then, the system time-varying state-space equation is established as

δ ⎤⎡ ⎤ y (t ) = Cx (t ) = ⎡ H H ⎦ ⎣ δ˙ ⎦ ⎣

(20)

Ξk(M) −i

and parameter p is required to satisfy p ≥ n/ m . Then, Eq. (20) is further rewritten in the following matrix form:

⎣

x˙ (t ) = Ac (t ) x (t ) + Bc (t ) u (t )

∑ Ξk(M) −i z k−i i=1

where

⎡ T F (t ) M (t ) = ⎢ 1 ⎢ T⋮ ⎣ FN (t )

(17)

m × 1

and ek ∈ denote the Kalman gain matrix and where Kk ∈ white innovation sequence, respectively. In addition, Ak ∈ n × n , Bk ∈ n × r , and Ck ∈ m × n are the discretized system, input, and output matrices, respectively. Then, Eqs. (16) and (17) are written as

(9)

i

(16)

The classical subspace methods for estimating the state vector x k are typically based on SVD or QR decomposition, and the original RPBSID method is required to construct a generalized Hankel matrix at each time instant to determine the state vector. However, all these methods require high computational complexity. To reduce computational cost, the state vector xk is estimated recursively in the current study. For the system in Eqs. (18)-(19) without the measurement noise ek , the output signals {yk − p , yk − p + 1 , ⋯, yk } from time k − p to k are written as

(15)

and Φ is the corresponding modal matrix. In numerical simulations, the modal matrices are obtained via finite element modeling. In actual onorbit identification, the system output signals can be measured using the gyroscopes on the center rigid body and the sensors on the appendages.

yk − p = Ck − p xk − p

(25)

yk − p + 1 = Ck − p + 1 A˜ k − p xk − p + Ck − p + 1 B˜k − p z k − p

(26)

yk = Ck A˜ k − 1 …A˜ k − p xk − p + Ck A˜ k − 1 …B˜k − p z k − p + ⋯+Ck B˜k − 1 z k − 1

(27)

3. Improved RPBSID algorithm In this section, the discrete innovation forms of Eqs. (13)-(14) are constructed and rewritten using a vector autoregressive with exogenous (VARX) predictor. An improved RPBSID algorithm is developed based on APSA, and the time-varying state-space model parameters are determined from the system's IeO data sequences.

Then, Eqs. (25)-(27) are expressed in the following matrix form:

y¯k − p = Γk − p xk − p + Hk − p z¯k − p

(28)

where

3.1. Recursive estimation of VARX parameter matrix

T T y¯k − p = ⎡ykT− p , ykT− p + 1 , ⋯, ykT ⎤ , z¯k − p = [z kT− p, z kT− p + 1, ⋯, z kT] ⎢ ⎥ ⎣ ⎦

The discretized innovation forms with noise for Eqs. (13) and (14) are expressed as follows: 159

Acta Astronautica 163 (2019) 157–167

Z. Ni et al.

Γk − p

similar purpose as the step-size parameter in conventional adaptive algorithms. Hence, we replace δd with a step-size parameter, i.e., μ , following the conventions. The general recursive form of the state vector xk at each time instant is obtained as follows:

Ck − p 0 ⎤ ⎡ ⎡ ⎤ ~ ~ ⎥ ⎢ Ck − p + 1 A ⎢ Ck − p + 1 Bk − p ⋱ ⎥ k−p ⎥ ⎢ , Hk − p = ⎢ = ⎥ ⋮ ⋯ 0 ⋮ ⎥ ⎢ ⎢ ~ ⎥ ~ ~ ~ ~ ⎥ ⎢C A ⎢ … ⋯ Ck Bk − 1 0 ⎥ …A ⎣Ck Ak − 1 Bk − p ⎦ ⎣ k k−1 k−p ⎦

xk = xk − 1 + μ

Notably, the elements of each block row in matrix Hk − p ∈ (p + 1) m × (p + 1)(m + r ) can be obtained through the recursive computation of the VARX parameter matrix Ξk in Eq. (23). Therefore, the impulse response matrix Hk−p is considered as known. If we define τk − p = y¯k − p − Hk − p z¯k − p , then Eq. (28) at the k − p time instant can be rewritten in the following least squares form:

τk − p = Γk − p xk − p

(41) (42)

We can choose 0 < μ ≪ 1 to guarantee the stability of the algorithm and to achieve a small steady-state misalignment because step-size parameter μ is derived from the minimum disturbance constraint δd , which should be sufficiently smaller than one to ensure computation convergence. Moreover, ν represents the regularization parameter, which is a selected positive number. In the recursive iteration of Eq. (41), the initial value Γ0 of matrix Γk − 1 should be a full-rank matrix. If no prior knowledge exists, then matrix Γ0 can be generally selected as follows:

APSA is then implemented now. For the k − 1 time instant, Eq. (29) is expressed in the following least squares form: (30)

where

τk − 1 = y¯k − 1 − Hk − 1 z¯k − 1

+ν

Σk − 1 = τk − 1 − Γk − 1 xk − 1

(29)

τk − 1 = Γk − 1 xk − 1

ΓkT− 1 sgn(Σk − 1) sgn(ΣkT− 1) Γk − 1 ΓkT− 1 sgn(Σk − 1)

(31)

I Γ0 = ⎡ n × n ⎤ ⎢ 0 ⎥ ⎦ ⎣

For Eq. (30), if we define a priori and a posteriori error vectors Σk−1 and Σk(p−)1, respectively, then for time instant k − 1,

(43)

Σk−1 = τk − 1 − Γk − 1 xk − 1

(32)

Similar to Eqs. 23 and 24, matrix Γk in Eq. (41) for each time instant is determined using the recursive least squares method as follows:

Σk(p−)1 = τk − 1 − Γk − 1 xk

(33)

Wk =

xkT Lk − 1 β2 + xkT Lk − 1 xk

(44)

Lk =

1 (Lk − 1 − Lk − 1 xk Wk ) β2

(45)

Subsequently, by minimizing the L1-norm of the posteriori error vector Σk(p−)1, the following minimum disturbance constraint problem can be obtained:

min ‖τk − 1 − Γk − 1 xk ‖1

Γk = Γk − 1 + (τk − Γk − 1 xk ) Wk

(34)

xk

subject to (‖xk − xk − 1 ‖2 )2 ≤ δd2

where β2 is a forgetting factor, and 0 ≪ β2 ≤ 1. The initial value L0 is selected as a unit matrix. Thus, the state vector xk for each time instant can be computed recursively using Eqs. (41)-(46).

(35)

where δd is a minimum disturbance parameter that ensures that the updating weight coefficient vector will not change dramatically. Then, using the Lagrange multiplier method, the unconstrained cost function can be obtained by combining Eqs. (34) and (35) as follows:

J (xk ) = ‖Σk(p−)1‖1 + βL (‖xk − x x − 1 ‖22 − δd2)

3.3. Recursive identification of time-varying state-space model parameters After the state vector xk has been determined using the system's IeO data, Eqs. (16) and (17) are rewritten as follows:

(36)

where βL is a Lagrange multiplier. The derivative of the cost function (36) with respect to the state vector xk is

∂J (xk ) = −ΓkT− 1 sgn(Σk(p−)1) + 2βL (xk − xk − 1) ∂xk

1 T Γk − 1 sgn(Σk(p−)1) 2βL

(37)

The corresponding system model parameter matrices Θk(x ) and Θk(y) are described as

(38)

Θk(x ) = [A˜ k , B˜k ]

(49)

Θk(y) = Ck

(50)

Using vector εk and output data yk , the recursive form of matrix Θk(y) is obtained via recursive least squares filtering as follows [31,38]:

(39)

Δk =

δd ΓkT− 1 sgn(Σk(p−)1) (p) T sgn(Σk(p−)T 1 ) Γk − 1 Γk − 1 sgn(Σk − 1)

(48)

T

When Eq. (39) is substituted into Eq. (38), the update equation for the state vector xk is

xk = xk − 1 +

yk = Θk(y) εk + ek

γk = [xkT , z kT] , εk = xk

δd (p) T sgn(Σk(p−)T 1 ) Γk − 1 Γk − 1 sgn(Σk − 1)

(47)

k

where

When Eq. (38) is substituted into Eq. (35), we obtain

1 = 2βL

= Θk(x ) γk

x

and the notation ‘sgn’ denotes a vector that is composed of the signs of the target vector. If we let the derivative of J (xk ) be equal to zero, then the following relationship can be derived:

xk = xk − 1 +

(46)

Θk(y) (40)

1 1 −1 Δk − 1 − Δk − 1 εk − 1 (β3 + εkT− 1 Δk − 1 εk − 1) (εkT− 1 Δk − 1 ) β3 β3

=

Θk(y−)1

+ (yk − 1 −

Θk(y−)1εk − 1)(εkT− 1 Δk )

(51) (52)

where the initial value Δ0 is typically selected as Δ0 = (1/ α2 ) I , α2 > 0 , and 0 ≪ β3 ≤ 1. The innovation noise sequence ek in Eq. (48) is then computed by

In Eq. (33), we use the priori error vector Σk − 1 to approximate Σk(p−)1 because the posteriori error vector Σk(p−)1 depends on the state vector xk , which is inaccessible before the current update. Furthermore, the minimum disturbance δd in Eq. (40) controls the convergence level of the algorithm, and thus should be considerably smaller than one to guarantee convergence. That is, the minimum disturbance δd serves a

ek = yk − Θk(y) εk matrix Θk(y)

(53)

is determined using Eq. (52), Eq. (53) is applied to When T update vector ek . Therefore, state vector γk = [xkT , z kT] can also be 160

Acta Astronautica 163 (2019) 157–167

Z. Ni et al.

updated. Similar to the computation of Θk(y) , matrix Θk(x ) is obtained recursively as follows:

Πk =

1 1 −1 Πk − 1 − Πk − 1 γk − 1 (β4 + γkT− 1 Πk − 1 γk − 1) (γkT− 1 Πk − 1) β4 β4

matrices and the singular value decomposition. Therefore, the computational cost of the improved method is lower than those of the classical RPBSID and SVD-based approaches. In particular, if the order n of the spacecraft system is high, then the advantage of the computational efficiency of the improved algorithm is even more apparent.

(54)

Θk(x ) = Θk(x−)1 + (xk − Θk(x−)1γk − 1)(γkT− 1 Πk )

(55) 4. Numerical simulations

where the initial value Π0 should satisfy Π0 = (1/ α3 ) I , α3 > 0 , and 0 ≪ β4 ≤ 1. After matrices Θk(x ) and Θk(y) are determined at each time instant, A˜ k , B˜k , and Ck are derived from Eqs. (49) and (50) by extracting the corresponding block matrices of Θk(x ) and Θk(y) . Matrices Ak , Bk , and Kk can be computed given that A˜ k = Ak − Kk Ck and B˜k = [Bk , Kk ]. 3.4. Computational complexity of the proposed method

The numerical model of the Engineering Test Satellite-VIII (ETSVIII) is established in this section. The time-varying state-space model parameters of the spacecraft are identified using the designed system's IeO signals. Moreover, the identification accuracy and computational efficiency of the proposed method, classical RPBSID algorithm, and SVD-based approach are compared.

The recursive procedures of the proposed algorithm can be summarized as follows

4.1. Parameters of the ETS-VIII model

(a) The VARX parameter matrix Ξk is constructed and updated using Eqs. (23) and (24). (b) Each block matrix of the impulse response matrix Hk − 1 in Eq. (31) is calculated via the recursive computation of the VARX parameter matrix Ξk . (c) The state vector xk is estimated recursively using Eqs. (41)-(46). (d) Parameter matrices Θk(x ) and Θk(y) are computed recursively using Eqs. (51)-(55). (e) The time-varying model matrices A˜ k , B˜k , and Ck are determined by extracting the corresponding block matrices of Θk(x ) and Θk(y) in Eqs. (49)-(50). (f) The state-space model matrices Ak , Bk , and Kk are computed from matrices A˜ k and B˜k .

The ETS-VIII was launched by Japan in 2006 to provide digital communications for mobile telephones and other mobile devices. The satellite has four large flexible appendages, namely, a pair of deployable antenna reflectors and a pair of solar panels. The satellite model is simplified based on the following conditions. The antenna reflectors are considered to be a plane truss structure, and the center rigid body of the satellite is regarded as a solid cuboid. The panels and reflectors are hinged with the center body by a link, and the four appendages are composed of a homogeneous material. The gravitational gradient torque is disregarded because ETS-VIII is a geosynchronous satellite. The entire satellite's center of mass is selected as the coordinate origin, and the origin of each appendage coordinate is established at the hinge joint that connects the center body and an appendage. The simplified configuration of ETS-VIII is shown in Fig. 1, where notations {s1, s2, a1, a2} denote the north/south solar panel and the A/B antenna reflector, respectively. The satellite model parameters used in the simulation are provided in Appendix A. The first eight orders of the appendages' frequencies obtained via finite element analysis are presented in Table 2, and the damping ratios of the appendages are ζs1 = ζs2 = ζa1 = ζa2 = 0.01. To simplify the discussion, if we select only the first two frequencies for each appendage in the model, then the state-space model order of the satellite system is n =(3 + 2 × 4) × 2 = 22.

The computational complexity of the improved RPBSID algorithm in identifying the time-varying model parameters is summarized in Table 1. We consider p = n/ m because the algorithm should satisfy p ≥ n/ m . If the system output dimension m = 1, then n = p ; otherwise, for m > 1, it has n > p. Therefore, we can use system order n instead of parameter p to estimate the computational complexity of the algorithm in Table 1. The computational complexity of the algorithm can be approximated as O (n2) flops per iteration for an n-order state-space system. The computational costs of the classical RPBSID and SVD-based methods should be at least O (n3) due to the multiplication of multiple

Table 1 Computational procedure of the improved RPBSID method, and the corresponding computational complexities. T

Initial conditions: Z0 = (1/ α1) Ip (m + r ) , Ξ0 = [ Im 0m × r … Im 0m × r ] , Γ0 = ⎡ InT 0nT× [n − m (p + 1)] ⎤ , L0 = In × n , x 0 = In × 1, Δ0 = (1/ α2 ) In , Π0 = (1/ α3 ) In + r + m , Θ0(x ) = [In , In × r , In × m], ⎢ ⎥    ⎣ ⎦ p (m + r )

Θ0(y ) = Im × n , α1,2,3 > 0 , 0 < μ ≪ 1, 0 < < β1,2,3,4 ≤ 1, and p ≥ n/ m . ν is a positive number, and Ηk − 1 in Step (3) is constructed from the VARX matrix Ξ . Recursive update: when k = 1, 2, … (1) Zk = (Zk − 1 − Zk − 1 ϕk (β1 + φkT Zk − 1 φk )−1 (φkT Zk − 1))/ β1

Computational complexity n (m + r )(2 + 5mn + 5rn)

(2)

Ξk = Ξk − 1 + (yk − Ξk − 1 φk )(φkT Zk )

n (m + r )(2m + mn + rn)

(3)

τ k − 1 = y¯k − 1 − Hk − 1 z¯k − 1

O (n2)

(4)

Σk − 1 = τ k − 1 − Γ k − 1 x k − 1

mn2 + mn

(5)

3mn2 + 4mn + n + m

(6)

x k = x k − 1 + μΓkT− 1 sgn(Σk − 1)/ sgn(ΣkT− 1) Γ k − 1 ΓkT− 1 sgn(Σk − 1) Wk = xkT Lk − 1/(β2 + xkT Lk − 1 x k )

(7)

Lk = (Lk − 1 − Lk − 1 x k Wk )/ β2

(8)

Γ k − 1 = Γ k − 1 + (τ k − Γ k − 1 x k ) Wk

+ν

2n2 + 2n

3n2 2mn2 + 2mn

(9)

Δk = (Δk − 1 − Δk − 1 εk − 1 (β3 +

(10)

Θk(y ) = Θk(y−)1 + (yk − 1 − Θk(y−)1εk − 1)(εkT− 1 Δk )

n2 + 2mn

(11)

ek = yk − Θk(y ) εk

mn

−1 εkT− 1 Δk − 1 εk − 1) (εkT− 1 Δk − 1))/ β3

(12)

Πk = (Πk − 1 − Πk − 1 γk − 1 (β4 + γkT− 1 Πk − 1 γk − 1)

(13)

Θk(x ) = Θk(x−)1 + (x k − Θk(x−)1γk − 1)(γkT− 1 Πk )

−1

(γkT− 1 Πk − 1))/ β4

5n2 + 2n

O (n2) 3n2 + 4(m + r ) n + (m + r )2

161

Acta Astronautica 163 (2019) 157–167

Z. Ni et al.

Fig. 3. Placement of sensors on each appendage. Fig. 1. Simplified ETS-VIII configuration. Table 2 First eight orders of frequencies for solar panels s1/s2 and antenna reflectors a1/a2 (unit: Hz). System order

Frequency of panels

Frequency of antennas

1 2 3 4 5 6 7 8

0.0767 0.5063 0.9958 1.4848 3.0459 3.0609 5.2750 5.3488

0.0887 0.2364 0.3698 0.4935 0.9236 0.9936 1.3815 1.7823

4.2. Design of IeO signals of the satellite Fig. 4. Spacecraft structure docking.

To simulate the control torque signal produced by the satellite reaction wheel, the designed input signals u = [u x , u y , uz ]T are shown in Fig. 2. In addition, the out-of-plane vibration response signals of certain nodes on each appendage are selected as the output vibration signals (Fig. 3).

4.3.1. Case 1: docking problem The following docking problem is considered in Case 1. The center rigid body of the spacecraft docks with two additional cube structures at the Z-axis (namely, the Yaw axis; Fig. 4). To simplify the modeling and analysis procedures in the simulation, we directly assume that the docking will not change the mass center of the satellite. Nevertheless, the position of the mass center will change in a typical docking operation because a case in which two space structures are docked symmetrically is rare. Therefore, the situation in this simulation is simply ideal, and the sole purpose is to verify the identification ability of the proposed recursive algorithm. The dimensions and mass of the docked cube structures are 2.35 m × 2.45 m × 3 m and 600 kg, respectively. Consequently, the spacecraft central mass mr changes with time t during the docking process as follows:

4.3. Results of identifying the parameters of the state-space model In the simulation, the measurement noise is a stationary zero-mean Gaussian random noise, and the selected signal-to-noise ratio is 40 dB. The parameters of the recursive method are given as follows: system sampling interval Δt = 0.1 s and Hankel matrix parameter p = 15 to ensure that the rank of the stacked matrix is larger than the system order n. In addition, the forgetting factors are β1 = β2 = β3 = β4 = 0.98, and the parameters are α1 = α2 = α3 = 0.9, μ = 0.01, and ν = 1. Thus, the following two cases of the time-varying model parameters are discussed.

Fig. 2. Designed satellite input signals. 162

Acta Astronautica 163 (2019) 157–167

Z. Ni et al.

Fig. 5. Designed test input signals: (a) torque input u x , (b) torque input u y , and (c) torque input uz .

mr =

⎧ 2400 kg 3000 kg ⎨ ⎩3600 kg

0≤ t <100 s 100 ≤ t <180 s 180 ≤ t ≤ 250 s

relative errors of the attitude angle responses of the three approaches are shown in Table 3. The results indicate that the three approaches have nearly the same computational accuracy. The relative errors of the proposed algorithm are only slightly larger than those of the other two methods. Further verification of additional simulation results proves that the identification accuracy of the improved method is basically the same as that of the original algorithm, and the relative errors of the three methods are less than 5% under this simulation condition. The relative error results prove that the proposed method can effectively identify the state-space model of the spacecraft system.

(56)

The docking process causes the dimensions and mass of the spacecraft's central rigid body to change. Therefore, the system state-space model parameters are also time-varying. Subsequently, the state-space model parameters of the satellite can be obtained using the improved RPBSID algorithm based on the IeO signals designed in Section 4.2. In addition to the proposed approach, the classical RPBSID and SVD-based methods are applied to this simulation [23,31]. Notably, a system has infinite sets of state-space model parameters {Ai , Bi , Ci} (i = 1,2, ..., ∞) and their specific element values vary. However, the different system state-space models should satisfy the same IeO relationship; thus, we can verify the accuracy of the identified time-varying model parameters by comparing the system response results. In particular, the test inputs {u x = sin(0.1t ) , u y = sin(0.05t ) , and uz = −sin(0.08t ) } (Fig. 5) are applied to the original and identified state-space models with zero initial state condition. The system response values that correspond to the attitude angles ψ , angular velocities ω , and the first three modal displacements η are then compared in Figs. 6 and 7. Although certain deviations occur in the final stage (Fig. 6), the overall test response values of the estimated state-space model are essentially consistent with those of the original model. The test response results of the estimated state-space model are generally consistent with those of the original model, and the results in Figs. 6 and 7 illustrate that all three methods can adequately identify the satellite's time-varying model parameters. Furthermore, the average

4.3.2. Case 2: rotation of the solar panels The solar panels of an on-orbit satellite always rotate around the pitch axis, such that they continually face the sun. Previous research results have demonstrated that the solar panel rotation of a satellite can cause a maximum 25% change in structure parameter [39]; thus, a satellite model is considered an LTV system. In Case 2, the solar panels are rotated at a constant angular velocity θ˙ , as shown in Fig. 8, where the rotation speed θ˙ = 0.2 deg/s in the numerical simulation. The other simulation conditions are the same as those in Case 1. Then, using the three aforementioned methods, the test response results of the identified state-space model are shown in Figs. 9 and 10. The corresponding average relative errors are provided in Table 4. Figs. 9 and 10 indicate that the attitude angle and angular velocity responses of the estimated state-space model are generally consistent with those of the original model and that the computed modal displacement results exhibit a certain delay. Table 4 shows that the relative error results of most variables of the improved method are less than those of the original algorithm. Case 2 also proves that the

Fig. 6. Comparison of the satellite attitude angles and angular velocities computed using the original and identified state-space models. 163

Acta Astronautica 163 (2019) 157–167

Z. Ni et al.

Fig. 7. Comparison of the satellite's first three modal displacements η computed by the original and identified state-space models.

computational accuracy of the two methods is nearly the same and can be acceptable. Although the relative errors of the two recursive methods may still be higher than that of the SVD-based approach in a noisy environment, certain signal de-noising technologies can be used to decrease the effect of noises and improve identification precision in future works.

Table 3 Average relative errors of the test output identified using the three methods. Attitude angle

Classical RPBSID (%)

SVD-based method (%)

Improved RPBSID (%)

ψx ψy

3.0835 2.3342

2.6808 2.4893

3.4380 3.2092

ψz

2.7766

3.6031

2.9923

4.4. Comparison of test response between the completely rigid and rigidflexible coupling spacecraft models In the previous section, the three methods are used to identify the state-space model of the satellite, and the corresponding test response results are compared. For the established rigid-flexible coupling spacecraft model in Section 2, if the flexible vibration influence of the spacecraft appendages is disregarded in the identification, that is, the obtained system output signals in Section 4.2 only include the spacecraft attitude signal but not the vibration signals of the appendages, the spacecraft is regarded as a complete rigid-body structure to identify the relevant time-varying state-space model. Then, by using the improved RPBSID method, the test response results of the completely rigid-body model are shown in Fig. 11. The test response results of the attitude angles and angular velocities in Fig. 11 validate that if a rigid-flexible coupling spacecraft system is directly considered as a rigid-body model to be identified without involving the flexible vibration effect, it exhibits an apparent distinction between the original and obtained state-space model parameters. Therefore, for certain spacecraft with obvious rigid-flexible coupling characteristics, to accurately identify the system model parameter, it is very important to obtain the structural vibration signals in

Fig. 8. Rotation of the solar panels.

Fig. 9. Comparison of satellite attitude angles and angular velocities computed using the original and identified state-space models. 164

Acta Astronautica 163 (2019) 157–167

Z. Ni et al.

Fig. 10. Comparison of the satellite's first three modal displacements η computed using the original and identified state-space models.

the computational efficiency of the three methods is further studied using the simulation results. Table 5 presents the computation time for simulations implemented in MATLAB for different system orders (30 Monte-Carlo experiments for each method). The results in Table 5 show that the computational efficiency of the classical and improved RPBSID algorithms is higher than that of the SVD-based method. The simulation results also validate the conclusion of the computational complexity analysis in Section 3.4; that is, as system order increases, the differences in computational time for these methods become increasingly apparent.

Table 4 Average relative errors of the test output identified using the three methods. Classical RPBSID (%)

SVD-based method (%)

Improved RPBSID (%)

Attitude angle ψx Attitude angle ψy

3.4729 1.6783

1.9646 1.2601

2.9664 1.6217

Attitude angle ψz

2.9235

2.6480

3.5811

Angular velocity ωx Angular velocity ω y Angular velocity ωz

5.8796 6.6838 7.4765

5.0194 6.8414 8.6877

5.7382 8.9119 7.7325

5. Conclusions To obtain the latest spacecraft model parameters in real time, an improved RPBSID method is proposed in this study and applied to identify the time-varying state-space model parameters of a spacecraft. The system state vector is computed using APSA and the recursive least squares method, thereby avoiding the requirement for a large amount of data, which is a serious issue in the original method. Accordingly, the time-varying state-space model is obtained recursively. The identification efficiency of the proposed method is proven theoretically by analyzing computational complexity. After the system state-space models are determined, the same test inputs are applied to the original and identified model parameters, and the system response values are compared to verify the validity of the improved recursive method. Two cases are studied via numerical simulations, in which the spacecraft model parameters may change. The response results show that the improved RPBSID method can effectively identify the timevarying state-space model matrices and have satisfactory identification accuracy. Further verification of additional simulation results proves

the identification process.

4.5. Comparison of computational efficiency for different system orders Finally, the identification accuracy and computational efficiency of the improved RPBSID method for different system model orders are examined. Three situations, in which the system order n of the established satellite model is set as 30, 38 and 46, are considered, with 30 Monte-Carlo experiments implemented for each situation. The comparison results of attitude angle ψx , angular velocity ωx , 1st modal displacement, and 1st modal velocity are shown in Fig. 12(a) and (b), (c), and (d), respectively. Fig. 12 indicates that the state-space model parameters can still be determined accurately when system order increases. The computational complexities of the proposed method, classical RPBSID, and SVD-based approaches are discussed in Section 3.4. Here,

Fig. 11. Comparison of the attitude angles and angular velocities for the rigid and rigid-flexible coupling spacecraft models using the improved RPBSID method. 165

Acta Astronautica 163 (2019) 157–167

Z. Ni et al.

Fig. 12. Comparison of the attitude angle ψx , angular velocity ωx , 1st modal displacement, and 1st modal velocity for the original and identified state-space models using the improved RPBSID method.

to identify the system model parameter accurately. The rigid-flexible coupling characteristics of large spacecraft in the modeling and identification procedures should not be disregarded. Finally, the comparison results of the computation times illustrate that the proposed algorithm achieve higher computational efficiency than the classical RPBSID and SVD-based methods. In particular, the simulation results validate the conclusion of the computational complexity theory analysis: when system order is high, the superior computational efficiency of the improved method is apparent. Nevertheless, similar to the classical RPBSID algorithm, the improved method requires prior knowledge, such as system model order.

Table 5 Average computation time in MATLAB of the three methods for different system orders n (unit: s). System order

Classical RPBSID

Improved RPBSID

SVD-based method

n = 22 n = 30 n = 38 n = 46

14.64 23.35 34.23 46.96

12.27 19.56 28.99 39.22

26.15 37.01 50.51 68.10

that the identification accuracy of the improved method is basically the same as that of the original algorithm, and the relative errors of the three methods are less than 5% under this simulation condition. Moreover, the spacecraft is regarded as complete rigid-body or rigidflexible coupling systems, respectively. The test response results of the attitude angles and angular velocities validate that for certain spacecraft with significant rigid-flexible coupling characteristics, obtaining the structural vibration signals in the identification process is important

Acknowledgment This work was supported by the National Natural Science Foundation of China (11502040, 51775541), the Key Program of the Chinese Academy of Sciences (Y4A3210301) and the Postdoctoral Science Foundation of China (2016M601354).

Appendix A. Simulation parameters of the ETS-VIII spacecraft The derivation of the time-varying matrix J (t ) in Case 2 is described by:

J (t ) = Jr − mr r˜r r˜r +

∑ TiT (t ) Ji Ti (t )− ∑ mi r˜i r˜i, i

i = {s1, s2, a1, a2} (A.1)

i

where subscript r denotes the satellite's central rigid body. The other simulation parameters are as follows: mr = 2400 kg, ms1 = ms2 = 115.2 kg, ma1 = ma2 = 165 kg

0 853.8 ⎤ 0 ⎤ ⎡ 3611.7 0 ⎡1683.4 Js1 = Js2 = ⎢ 0 3226.8 0 ⎥ 18.2 0 ⎥, Ja1 = Ja2 = ⎢ 0 0 3291.3⎦ 0 3629.9 ⎦ ⎣ 0 ⎣ 853.8 0 0 ⎤ − 0.4827 0 ⎤ ⎡11858.5 ⎡ 0 Jr = ⎢ 0 11762.5 0 ⎥, r˜r = ⎢ 0.4827 0 0 ⎥, r˜s1 = 0 2305⎦ 0 0⎦ ⎣ 0 ⎣ 0

− 0.4827 5 ⎤ ⎡ 0 r˜s2 = ⎢ 0.4827 0 0 ⎥, r˜a1 = 0 0⎦ ⎣ −5

− 0.4827 − 5⎤ ⎡ 0 0 0 ⎥ ⎢ 0.4827 0 0 ⎦ ⎣ 5

0 ⎤ − 0.4827 0 ⎤ − 0.4827 ⎡ 0 ⎡ 0 0 1.175⎥ 0 − 1.175⎥, r˜a2 = ⎢ 0.4827 ⎢ 0.4827 0 ⎦ 1.175 0 ⎦ − 1.175 ⎣ 0 ⎣ 0

˙ ) 0 sin(θt ˙ )⎤ ˙ ) 0 − sin(θt ˙ )⎤ ⎡− cos(θt ⎡ cos(θt 1 0 0⎤ ⎥, Ta1 = Ta2 = ⎡ Ts1 = ⎢ 0 −1 0 ⎥, Ts2 = ⎢ 0 1 0 0 1 0⎥ ⎢ ⎢ ⎢ 0 0 1⎦ ˙ ) ˙ )⎥ ˙ ) 0 cos(θt ˙ ) ⎥ ⎣ θt θt θt sin( 0 cos( sin( ⎣ ⎦ ⎣ ⎦ 166

Acta Astronautica 163 (2019) 157–167

Z. Ni et al.

[20] I. Yamaguchi, T. Kida, T. Kasai, Experimental demonstration of LSS system identification by eigensystem realization algorithm, Proceedings of the 1995 American Control Conference. Seattle, WA, USA, June 21–23, 1995. [21] I. Yamaguchi, T. Kida, System identification experiments of large space structures, Trans. Soc. Instrum. Control Eng. E-1 (1) (2001) 312–321. [22] Y. Xie, P. Liu, G.P. Cai, On-orbit identification of flexible parameters of spacecraft, Proc. Inst. Mech. Eng. - Part K J. Multi-body Dyn. 230 (2) (2016) 191–206. [23] K.F. Liu, L. Deng, Identification of pseudo-natural frequencies of an axially moving cantilever beam using a subspace-based algorithm, Mech. Syst. Signal Process. 20 (1) (2006) 94–113. [24] M. Majji, J.N. Juang, J.L. Junkins, Time-varying eigensystem realization algorithm, J. Guid. Contr. Dynam. 33 (1) (2010) 13–28. [25] I. Markovsky, J. Goos, K. Usevich, et al., Realization and identification of autonomous linear periodically time-varying systems, Automatica 50 (6) (2014) 1632–1640. [26] A. Estes, M. Majji, J.N. Juang, Time-varying methods for identification of constrained flexible structures, AIAA/AAS Astrodynamics Specialist Conference. San Diego, CA, USA, August 4–7, 2014. [27] K.P. Yu, K. Yang, Y. Bai, Experimental investigation on the time-varying modal parameters of a trapezoidal plate in temperature-varying environments by subspace tracking-based method, J. Vib. Contr. 21 (16) (2015) 3305–3319. [28] M. Arjomandi-Lari, M. Karimi, Generalized YAST algorithm for signal subspace tracking, Signal Process. 117 (2015) 82–95. [29] F. Haddadi, Steady-state statistical performance analysis of subspace tracking methods, IEEE Trans. Signal Process. 64 (18) (2016) 4781–4791. [30] Z.Y. Ni, R.N. Mu, G.B. Xun, et al., Time-varying modal parameters identification of a spacecraft with rotating flexible appendage by recursive algorithm, Acta Astronaut. 118 (2016) 49–61. [31] I. Houtzager, J.W. van Wingerden, M. Verhaegen, Recursive predictor-based subspace identification with application to the real-time closed-loop tracking of flutter, IEEE Trans. Contr. Syst. Technol. 20 (4) (2012) 934–949. [32] J.W. van Wingerden, A. Hulskamp, T. Barlas, et al., Two-degree-of-freedom active vibration control of a prototyped “smart” rotor, IEEE Trans. Contr. Syst. Technol. 19 (2) (2011) 284–296. [33] M. Bergamasco, M. Lovera, Continuous-time predictor-based subspace identification for helicopter dynamics, 37th European Rotorcraft Forum. Gallarate, Italy, September 13–15, 2011. [34] J. Wartmann, S. Seher-Weiss, Application of the predictor-based subspace identification method to rotorcraft system identification, 39th European Rotorcraft Forum. Moscow, Russia, September 3–6, 2013. [35] P.C. Hughes, Modal identities for elastic bodies, with application to vehicle dynamics and control, J. Appl. Mech. 47 (1) (1980) 177–184. [36] H.B. Hablani, Constrained and unconstrained modes-some modeling aspects of flexible spacecraft, J. Guid. Contr. Dynam. 5 (2) (1982) 164–173. [37] A.H. Sayed, Adaptive Filters, Wiley, New York, 2008. [38] S.R. Seydnejad, Real-time identification of time-varying ARMAX systems based on recursive update of its parameters, J. Dynamic Syst. Meas. Contr. Trans. ASME 136 (3) (2014) 696–706. [39] T. Nagashio, T. Kida, T. Ohtani, et al., Design and implementation of robust symmetric attitude controller for ETS-VIII spacecraft, Contr. Eng. Pract. 18 (12) (2010) 1440–1451.

References [1] J. Cooper, J. Wright, Spacecraft in-orbit identification using eigensystem realization methods, J. Guid. Contr. Dynam. 15 (2) (1992) 352–359. [2] I. Yamaguchi, T. Kida, K. Komatsu, et al., ETS-VI on-orbit system identification experiments, JSME Int. J. C Mech. Syst. Mach. Elem. Manuf. 40 (4) (1997) 623–629. [3] S. Adachi, I. Yamaguchi, T. Kida, et al., On-orbit system identification experiments on Engineering Test Satellite-VI, Contr. Eng. Pract. 7 (7) (1999) 831–841. [4] I. Yamaguchi, T. Kasai, H. Igawa, Multi-input multi-output system identification using impulse responses, J. Space Eng. 1 (1) (2008) 1–11. [5] T. Kasai, I. Yamaguchi, H. Igawa, et al., On-orbit system identification experiments of the Engineering Test Satellite-VIII, Trans. Jpn Soc. Aeronaut. Space Sci. Aerospace Technol. Jpn. 7 (26) (2009) 79–84. [6] T. Nagashio, T. Kida, Y. Hamada, et al., Robust two-degrees-of-freedom attitude controller design and flight test result for Engineering Test Satellite-VIII spacecraft, IEEE Trans. Contr. Syst. Technol. 22 (1) (2014) 157–168. [7] R.S. Pappa, J.N. Juang, Galileo spacecraft modal identification using an eigensystem realization algorithm, J. Astronaut. Sci. 33 (1) (1985) 15–33. [8] J.N. Juang, M. Phan, L.G. Horta, et al., Identification of observer/Kalman filter Markov parameters-theory and experiments, J. Guid. Contr. Dynam. 16 (2) (1993) 320–329. [9] P. Nuthi, K. Subbarao, Computational adaptive optimal control of spacecraft attitude dynamics with inertia matrix identification, J. Guid. Contr. Dynam. 40 (5) (2017) 1258–1262. [10] D. Mackison, Identification and adaptive control of a satellite flight control system, AIAA/AAS Astrodynamics Specialist Conference and Exhibit. Monterey, CA, USA, August 5–8, 2002. [11] S. Aoues, F.L. Cardoso-Ribeiro, D. Matignon, et al., Modeling and control of a rotating flexible spacecraft: a port-Hamiltonian approach, IEEE Trans. Contr. Syst. Technol. 99 (2017) 1–8. [12] Y. Xie, P. Liu, G.P. Cai, Modal parameter identification of flexible spacecraft using the covariance-driven stochastic subspace identification (SSI-COV) method, Acta Mech. Sin. 32 (4) (2016) 710–719. [13] D. Ivanov, S.V. Meus, A.V. Ovchinnikov, et al., Methods for the vibration determination and parameter identification of spacecraft with flexible structures, J. Comput. Syst. Sci. Int. 56 (2) (2017) 311–327. [14] H.K. Rad, H. Salarieh, A. Alasty, et al., Boundary control of anti-symmetric vibration of satellite with flexible appendages in planar motion with exponential stability, Acta Astronaut. 147 (2018) 219–230. [15] N. Ji, J.K. Liu, Vibration control for a flexible satellite with input constraint based on Nussbaum function via backstepping method, Aero. Sci. Technol. 77 (2018) 563–572. [16] A. Meguro, K. Shintate, M. Usui, et al., In-orbit deployment characteristics of large deployable antenna reflector onboard Engineering Test Satellite VIII, Acta Astronaut. 65 (9) (2009) 1306–1316. [17] D. Entekhabi, E.G. Njoku, P.E. O'Neill, et al., The soil moisture active passive (SMAP) mission, Proc. IEEE 98 (5) (2010) 704–716. [18] H.M. Kim, E.B. Bokhour, Mir structural dynamics experiment: a flight experiment development, The 38th AIAA Structures, Structural Dynamics, and Materials Conference. Kissimmee, FL, USA, April 7–10, 1997. [19] Y. Xie, P. Liu, G.P. Cai, On-orbit frequency identification of spacecraft based on attitude maneuver data, Aircraft Eng. Aero. Technol. 90 (1) (2018) 33–42.

167