State Trajectory Generation of MIMO Multirate
Feedforward for Perfect Tracking Control in
High-Precision Stage
Masahiro Mae and Hiroshi Fujimoto
The University of Tokyo
5-1-5, Kashiwanoha, Kashiwa, Chiba, 277-8561, Japan
Phone: +81-4-7136-3881
Email: mmae@ieee.org, fujimoto@k.u-tokyo.ac.jp

Abstract—Multirate feedforward control is widely used for thecontrol of precise systems such as the high-precision stages andthe hard disk drives. In the multirate feedforward control, it ispossible to obtain a control input that achieves perfect tracking
for the state trajectory. Therefore, it is necessary to generatea state trajectory corresponding to the output trajectory inadvance. In this paper, the authors propose a state trajectory generation from an output trajectory for a multi-input multi-output
system using singular value decomposition. The effectiveness of
this method is verified through simulation and experiments.
Index Terms—singular value decomposition, state trajectory
generation, multirate feedforward control, perfect tracking control, multi-input multi-output system, high-precision stage

I. I NTRODUCTION
Multirate feedforward control is widely used for the controlof the precise systems such as the high precision stages
and the hard disk drives. Multirate feedforward solves the
problem of discretize unstable zeros in the discrete system
and achieves perfect tracking control (PTC) [1]. The word
“perfect tracking control (PTC)” is defined as “the plant output
perfectly tracks the desired trajectory with zero tracking error
at every sampling point” [2].
The high-precision stages have an important role in the
manufacturing of semiconductors and liquid crystal panels[3]. Multirate feedforward control is commonly applied tothe high-precision stages. The high-performance demand of
semiconductors and liquid crystal panels is increasing year
by year in the world [4]. The conventional high-precision
stages are single-input single-output (SISO) 1 degree of freedom (DOF) systems mechanically constrained other than longstroke. However, in recent years, from the demands of further
precision, the high-precision stages become 6 DOF systems
which are floating by air bearings without mechanical frictions,and they are attempted to be controlled as multi-input multioutput (MIMO) systems [5]. In MIMO systems represented
by 6 DOF system, the coupling of each DOF is one ofthe problems for the high-precision control [6]. For coupled
MIMO systems, it is common to use a precompensator todecouple and then control by using SISO controllers, but theprecompensator is not considering the discretize unstable zeroproblems. Therefore, in order to achieve perfect tracking control, the multirate feedforward controller for MIMO systems
is needed to be designed. In multirate feedforward control for
MIMO systems, the design of the controller is not uniquedepending on the select the B matrix from the generalizable
controllability indices, and several kinds of controllers aredesigned [7]. Using multirate feedforward control for MIMO
systems, PTC can be achieved in several patterns [8].
In the multirate feedforward control, it is possible to obtain acontrol input for achieving the perfect tracking control for thestate trajectory. Therefore, it is necessary to calculate the state
trajectory corresponding to the output trajectory in advance.
In the conventional multirate feedforward control for the
SISO system, the state trajectory corresponding to the output
trajectory is calculated by convolving the output equation of
the controllable canonical form in the time domain. If thesystem has the unstable zeros in the continuous time domain,
it is not possible to directly calculate the stable inverse system.
In this case, the preactuation perfect tracking control (PPTC)
using time axis reversal [9] or the stable inversion [10] havebeen used. On the other hand, in the multirate feedforward
control for the MIMO system, it is difficult to calculate thestate trajectory corresponding to the output trajectory becauseof the presence of off-diagonal terms in the C matrix ofthe output equation, and the complexity of the controllable
canonical form in MIMO systems.
In this paper, the authors propose the method to generate
the state trajectory corresponding to the output trajectory for
multirate feedforward control in MIMO system using singularvalue decomposition (SVD). The authors control the highprecision stage with 6 DOF, and consider 2DOF system withthe translation along the x axis and the pitching around they axis. The authors propose the method of state trajectory
generation using singular value decomposition and the multirate feedforward controller design in the MIMO system, andthe effectiveness of the proposed method compared with theconventional method is verified through the simulation and theexperiment.

(b) Coupling problem between x and θy .

(a) Photograph of the 6-DOF high-precision
stage.

(c) Block diagram of the plant.

Fig. 1. Details of the plant.

Fig. 2. Block diagram of the controller with the plant which does not have unstable intrinsic zeros. S, H, and L denote a sampler, holder, and lifting operator[11], respectively, z and zs denote esTr and esTu , respectively.

II. C ONTROLLER D ESIGN
A. State Trajectory Generaion
In this paper, the authors think about the plant without thecontinuous unstable zeros.
1) State Trajectory Generation for Single-Input SingleOutput System: The conventional method of generating a state
trajectory from an output trajectory in the SISO system will
be described below.
The transfer function of the plant is given by
B(s)
A(s)
bm sm + bm−1 sm−1 + · · · + b0=
sn + an−1 sn−1 + · · · + a0

P (s) =

(1)

The state space representation of this plant in controllable
canonical form is given by 
x1 (t)0 0d  x2 (t)   . = dt  ..  
xn (t)−a0y(t)= b0

b1

1
0

···
···
..
.

−a1

···

···

bm

  
x1 (t)0  x2 (t)  0    ..  +  ..  u(t),  .  .
−an−1
xn (t)1
x1 (t)

 x2 (t) 
0 · · · 0  . .
.
 . 
0
0



(2)

(3)

xn (t)

Since each state variable is in a differential relation fromthe state equation (2), the convolution of the output equation(3) in the time domain gives the state trajectory xd (t) =T
x1d (t) x2d (t) · · · xnd (t)from the output trajectoryyd (t) as1−1
yd (s)(4)x1d (t) = L
bm sm−1 + · · · + b1 s + b0
where L−1 is the inverse Laplace transform.
2) State Trajectory Generation for Multi-Input MultiOutput System: The proposed method of generating the state
trajectory from the output trajectory in the MIMO system will
be described below.
This proposed method is calculated by the following 6 steps.
1) Laplace transform state equation and output equation
2) Using the singular value decomposition for B matrix
3) Divide the state equation into two parts (Wu (s) andWl (s)) with and without input u(s)4) Concatenate output equations and Wl (s) vertically5) Solve about state trajectory x(s)6) Calculate x(t) from inverse Laplace transform
The details of each step are explained below.
In this paper, the plant is assumed that m < n, rank(B) =m, and does not have the continuous unstable zeros.
The state equation and the output equation of m input m
output nth order plant are given byẋ(t) = Ax(t) + Bu(t),(5)y(t) = Cx(t),(6)

where the state variables x(t) ∈ Rn×1 , the input u(t) ∈Rm×1 , the output y(t) ∈ Rm×1 , and the matrices A ∈ Rn×n ,B ∈ Rn×m , and C ∈ Rm×n .
The singular value decomposition of the matrix B is given
byFrom (16), the state trajectory x̃(s) from the output trajectory y(s) in the MIMO system is given by

−1 

CU
y(s)
x̃(s) =.
(17)
Wl (s)O

B = U ΣV ∗ ,From the state variable conversion of (17) using the transformation matrix U as x(s) = U x̃(s), the state trajectoryx(s) from the output trajectory y(s) in the MIMO system isgiven by

−1 

CU
y(s)
(U x̃(s)) = U
Wl (s)
O

−1 

CU
y(s)
⇔ x(s) = U.
(18)
Wl (s)
O(7)

where U ∈ Rn×n , Σ ∈ Rn×m , V ∈ Rm×m .
The elements of Σ are given by
 
∆
Σ=
O
∆ = diag(σi ) (i = 1 · · · m ∈ N),
where σi (i = 1, 2, . . . , m ∈ N) are the singular values of thematrix B.
Multiply both sides of (7) by U −1 from the left
U −1 B = ΣV ∗


∆V ∗
=
.
O

(8)

Laplace transform of the state equation (5) is given bysx(s) = Ax(s) + Bu(s).

B. Multirate Feedforward(9)

The state variable conversion of (9) using the transformation
matrix U −1 as x̃(s) = U −1 x(s) is given bys(U −1 x(s)) = (U −1 AU )(U −1 x(s)) + (U −1 B)u(s)⇔ (sI − U −1 AU )x̃(s) = (U −1 B)u(s)∆V ∗⇔ W (s)x̃(s) =u(s).
OW (s) of (10) is represented as
Wu (s)
W (s) =,Wl (s)

By inverse Laplace transform of (18), the state trajectory
xd (t) from the output trajectory yd (t) in the MIMO system
is given by
" 
−1 
#
CU
yd (s)−1
xd (t) = L
U.
(19)
Wl (s)
O

(10)

(11)

1) Multirate Feedforward for Single-Input Single-Output
System: Multirate feedforward control can achieve PTC [1].
A digital tracking control system usually has two samplers
for the reference signal r(t) and the output y(t), and one
holder on the input u(t), as shown in Fig. 2. Therefore, there
exist three time periods Tr , Ty , and Tu which represent theperiods of r(t), y(t), and u(t), respectively, as shown in Fig.
3. Additionally, the longer period of Tr and Ty is defined as
the frame period Tf . In the case of the SISO nth order plant,
Tr = nTy = nTu = Tf .
Consider the continuous-time nth order plant described by
the state equation (20) and the output equation (21).
ẋ(t) = Ac x(t) + bc u(t)(20)y(t) = cc x(t)(21)

where W (s) = sI −U −1 AU , Wu (s) ∈ Cm×n , and Wl (s) ∈C(n−m)×n .
From (10) and (11)Wu (s)∆V ∗x̃(s) =
u(s).
(12)
O
Wl (s)

From discretization of (20) and (21) by zero order hold insampling period Tu , the discrete-time plant becomes the state
equation (22) and output equation (23).

The extracted lower (n − m) rows of (12) is given byx[k + 1] = As x[k] + bs u[k]Wl (s)x̃(s) = O.

(13)

Laplace transform of the output equation (6) is given byy(s) = Cx(s).

(14)

The state variable conversion of (14) using the transformation matrix U −1 as x̃(s) = U −1 x(s) is given byy(s) = (CU )(U −1 x(s))⇔ (CU )x̃(s) = y(s).

(15)

The concatenate matrix from (13) and (15) is given by




CU
y(s)
x̃(s) =.
(16)
Wl (s)
Oy[k] = cs x[k]

(22)
(23)

where the matrices As , bs , and cs are given by
As = eAc Tu ,
Z Tu
bs =
eAc τ bc dτ,

(24)
(25)

0

cs = cc .

(26)

By lifting the discrete-time state equation (22) and the output
equation (23), the state equation (27) and the output equation(28) are given byx[i + 1] = Ax[i] + Bu[i],y[i] = cx[i],

(27)
(28)

selected from{bc1 , · · · , bcm , Ac bc1 , · · · , Ac bcm , · · · , An−1bcm }.
c
Letting ϕ be a set of these n vectors, σl , and N are definedbymX

σl = number{k|Ak−1
bcl ∈ ϕ},c(40)

σl = n,(41)

N = max(σl ).

(42)

l=1

Fig. 3. Multirate sampling control at same interval.

where the matrices A, B, and c, and i are given by
A = Ans ,bs
B = An−1
s

(29)
An−2
bs
s

···

As bs



bs ,

(30)

c = cs ,(31)

x[i] = x(iTr ).

(32)(27) and (28) are given from the interval between t = iTr =kTu and t = (i + 1)Tr = (k + n)Tu . Input and output vectors
u[i] (33) and y[i] (34) are given byTu[i]= u1 [i] u2 [i] · · · un [i]T= u(kTu ) u((k + 1)Tu ) · · · u((k + n − 1)Tu ) ,Ty[i]= y1 [i] y2 [i] · · · yn [i]T= y(kTy ) y((k + 1)Ty ) · · · y((k + n − 1)Ty ) .

x[i + 1] = Ax[i] + Bu[i],uf f [i] = B−1(I − z−1(43)A)x[i + 1],(44)

where the matrices A, x[i], u[i], and z are given by
A = eAc Tf ,

(45)

(33)

x[i] = x(iTf ),

(46)

(34)


T
u[i]= u1 [i] · · · um [i]= u11 [i] · · · u1σ1 [i] u21 [i] · · ·From (27), the control input uf f [i] to achieve PTC are givenby
uf f [i] = B −1 (I − z −1 A)x[i + 1].

In the case of MIMO system, n (= plant order) number ofvectors are selected from generalized controllability indices,and the full row rank matrix B are designed. Therefore, thefeedforward controllers have different forms according to theirchoice.
From (43), the control input uf f [i] to achieve PTC are givenby (44).

ẋ(t) = Ac x(t) + Bc u(t)(36)y(t) = Cc x(t)Bc = bc1 · · · bcmT
Cc = cc1 · · · ccp

(37)
(38)
(39)

where the plant state x ∈ Rn , the plant input u ∈ Rm , andthe plant output u ∈ Rp .
Generalized controllability indices are defined as follows[7].
Definition (Generalized Controllability Indices). Generalized
controllability indices of (Ac , Bc ) are defined as follows forAc ∈ Rn×n and Bc = [bc1 , · · · , bcm ] ∈ Rn×m . If (Ac , Bc )is a controllable pair, n linearly independent vectors can be

T

,

(47)

,

(48)

Tf = N Tu .

(49)

z= e

(35)

The block diagram of the control system is shown in Fig.
2. L is a discrete-time lifting operator [11]. L−1 outputs theelements of nth dimensional vector uf f [i], which is input
every period Tr , in order from 1 to n by Tu = Tr /n.
2) Multirate Feedforward for Multi-Input Multi-Output System: In m-input p-output nth order MIMO system, the state
equation (36) and the output equation (37) of the continuoustime plant are given bysTf

umσm [i]III. E XPERIMENT
A. ModelingIn the simulation and experiment, the authors control thefine stage of the 6 DOF high-precision stage shown in Fig.
1(a). This fine stage is supported by 6 DOF air bearing gravity
canceller. In this paper, controlling 2 DOF of the translationx along the x axis and the pitching θy around and y axis, as
shown in Fig. 1(b).
The equations of motion of the translation and pitching of
the stage are given by (50) and (51) [5].
(Mx1 + Mx2 )ẍg1 + Cx1 ẋg1 + Kx1 xg1 + Mx2 Lg2 θ̈y = fx

(50)

(Mx2 L2g2 + Jθy )θ̈y + Cθy θ̇y + Kθy θy + Mx2 Lg2 (ẍg1 − gθy )= τy + fx Lf x

(51)

Convert xg1 to observable xm by (52).
xm (s) = xg1 (s) + Lm θy (s)(52)

From (50), (51), and (52), transfer functions g11 to g22 in Fig.
1(c) are given. In the proposed method, this plant model isrepresented in the balanced realization.
The parameters of the stage are shown in Table I, which is
given from fitting in the frequency domain, shown in Fig. 4.

Since the problem of the unstable discretization zero is notconsidered, PTC cannot theoretically be achieved.
In the proposed method, a dual-input dual-output multirate
feedforward controller is designed for a dual-input dual-output
plant as shown in Fig. 5(c).In this paper, the selection of
the generalized controllability indices in the proposed method(σ1 , σ2 ) = (2, 2) is given by
B = As bs1 bs1 As bs2 bs2(53)

TABLE IM ODEL PARAMETERS .
Symbol
xm
xg1
xg2
θy
fx
τy
Mx1
Cx1
Kx1
Mx2
Jθy
Cθy
Kθy
Lm
Lg2
Lf x

Meaning
Measured position of the fine stage
Position of the CoG of the planar air bearing and the air gyro
Position of the CoG of the fine stage
Measured attitude angle of the fine stage
Input force of the fine stage in the x direction
Input torque of the fine stage in the θy direction
Mass of the planar air bearing and the air gyro
Viscosity coefficient in the xg1 motion
Spring coefficient in the xg1 motion
Mass of the fine stage
Moment of inertia of the fine stage
Viscosity coefficient of the fine stage in the θy motion
Spring coefficient of the fine stage in the θy motion
Distance between the measurement point of xm and the CoRDistance between the CoR and the CoG of the fine stage
Distance between the CoR of the fine stage and the actuation point

-50

-50

-100

-100

-150

-150

-200
0
10

1

10

2

10

3

10

-200
0
10

180

180

90

90

0

0

-90

-90

-180
100

101

102

103

-180
100

(a) g11

Value
–
–
–
–
–
–
0.077 kg
200 N/(m/s)
5000 N/m
5.3 kg0.10 kgm2
1.6 Nm/(rad/s)
1200 Nm/rad−0.028 m−0.051 m−0.0026 m

10

1

10

2

10

3

101

102

103

0

-100

-50

-150

-100

C. Simulation Results
The simulation results are shown in Fig. 6. From the simulation results, PTC cannot be achieved with the conventional
method, but it can be achieved in the proposed method. The
effectiveness of the proposed method is verified.

(b) g12-50

In Fig. 4, G4 is used for MIMO plant model of the
proposed method, and G2 is used for SISO plant model ofthe conventional method 1 and 2.
In the simulation, the feedback controller Cf b [zs ] is 0 in
the all methods. In the experiment, Cf b [zs ] are PID controllers tuned to have 20 Hz closed roop pole for all 6 DOF(x, y, z, θx , θy , θz ) in the all methods.
7th order polynomial trajectory of 0 µm to 100 µm in 0 ms
to 20 ms is given for xref
m and 7th order polynomial trajectoryof 0 µrad to 100 µrad in 0 ms to 20 ms is given for θyref . Tu =200 µs, N = max(σ1 , σ2 ), and Tr = N Ty = N Tu = Tf .

D. Experiment Results
-200
0
10

1

10

2

10

3

10

-150
0
10

180

180

90

90

0

0

-90

-90

-180
100

101

(c) g21

102

103

-180
100

10

1

10

2

10

3

101

102

103

(d) g22

Fig. 4. Frequency responses of the plant. G4 and G2 denote a 4th order
MIMO nominal plant and a 2nd order SISO nominal plant, respectively.

B. Conditions
The block diagram of the simulation is shown in Fig. 2.
Generally, a dual-input dual-output plant with the interference
is represented by a block diagram as shown in Fig. 1(c). g11
to g22 represent the transfer function of each path.
Conventional method 1 ignores the interference by g12 andg21 and the multirate feedforward controllers are designed foreach SISO system g11 and g22 , respectively, as shown in Fig.
5(a). Coupling is suppressed by the feedback controller.
Conventional method 2 designs a precompensator forcontinuous-time plants that cancels the interference between
axes by g12 and g21 , and the multirate feedforward controllers
are designed for each decoupled SISO system g11 and g22 ,respectively, as shown in Fig. 5(b). In this method, theprecompensator is discretized by the bilinear transformation.

In the experiment, the proposed method of (σ1 , σ2 ) = (2, 2),
the conventional method 1, and the conventional method 2 arecompared under the same conditions as the simulation.
The experiment result is shown in Fig. 7. From the experiment results, the tracking errors of xm are close values in all
methods. The tracking errors of xm in the conventional method
2 and the proposed method become bigger than the simulation
results due to the effect of modeling error, shown in Fig. 4.
However, the tracking error of θy in the proposed method is
the smallest compared with those of the conventional methods.
The effectiveness of the proposed method is verified.
IV. C ONCLUSION
In this paper, the authors propose the multirate feedforward
controller and the generating method of the state trajectory
corresponding to the output trajectory for MIMO systems. In
the MIMO multirate feedforward controller design, the authors
show that several kinds of controllers can be designed for
the B matrix selected from the generalized controllability
indices. In addition, the authors show the generating methodof the state trajectory in any state space representation using
singular value decomposition. The effectiveness of the proposed method is verified by simulations and experiments intranslation and pitching motions of the 6 DOF high-precision
stage. The case in which the plant has continuous unstablezeros and the extension from 2 DOF to 6 DOF are futurework.

(a) Conventional method 1.

(b) Conventional method 2.

(c) Proposed method.

Fig. 5. Block diagram of the control methods. These show only feedforward controllers.
10-4

1

10-6

1

0.8

0

0.6

-1

0.4

-2

0.2

-3

10-4

3

10-5

0

2
-10

0

1

0

-4
0

0.005

0.01

0.015

0.02

0

0.005

(a) xm

0.01

0.015

0.02

-20

0

0.005

(b) exm

0.01

0.015

0.02

0

0.005

(c) θy

0.01

0.015

0.02

0.015

0.02

(d) eθyFig. 6. Simulation results.
10-5

10-6

10-5

20

0

10-5
5

10
15

8

-5

0

6

10
-10

4

-5

5
2
-15

0
0

0.005

0.01

0.015

0.02

(a) xm

0
0

0.005

0.01

0.015

0.02

-10
0

0.005

(b) exm

0.01(c) θy

0.015

0.02

0

0.005

0.01

(d) eθyFig. 7. Experimental results.

R EFERENCES[1] H. Fujimoto, Y. Hori, and A. Kawamura, “Perfect tracking control basedon multirate feedforward control with generalized sampling periods,”
IEEE Transactions on Industrial Electronics, vol. 48, no. 3, pp. 636–
644, jun 2001.
[2] M. Tomizuka, “Zero Phase Error Tracking Algorithm for Digital Control,” Journal of Dynamic Systems, Measurement, and Control, vol. 109,
no. 1, p. 65, 1987.
[3] T. Oomen, “Advanced Motion Control for Precision Mechatronics: Control, Identification, and Learning of Complex Systems,” IEEJ Journal ofIndustry Applications, vol. 7, no. 2, pp. 127–140, 2018.
[4] G. E. Moore, “Cramming More Components Onto Integrated Circuits,”
Electronics, vol. 38, no. 8, 1965.
[5] W. Ohnishi, H. Fujimoto, K. Sakata, K. Suzuki, and K. Saiki, “Decoupling Control Method for High-Precision Stages using Multiple
Actuators considering the Misalignment among the Actuation Point,
Center of Gravity, and Center of Rotation,” IEEJ Journal of Industry
Applications, vol. 5, no. 2, pp. 141–147, 2016.
[6] H. Butler, “Position Control in Lithographic Equipment,” IEEE Control
Systems, vol. 31, no. 5, pp. 28–47, oct 2011.
[7] H. Fujimoto, “General Framework of Multirate Sampling Control andApplications to Motion Control Systems,” Ph.D. dissertation, 2000.
[8] M. Mae, W. Ohnishi, H. Fujimoto, and Y. Hori, “Perfect Tracking
Control of Dual-Input Dual-Output System for High-Precision Stagein Translation and Pitching Motion,” Proceedings of the International
Workshop on Sensing, Actuation, Motion Control, and Optimization,
no. 5, pp. TT6–5, 2018.
[9] W. Ohnishi, T. Beauduin, and H. Fujimoto, “Preactuated multirate
feedforward for a high-precision stage with continuous time unstablezeros,” IFAC-PapersOnLine, vol. 50, no. 1, pp. 10 907–10 912, jul 2017.
[10] J. van Zundert, J. Bolder, S. Koekebakker, and T. Oomen, “Resourceefficient ILC for LTI/LTV systems through LQ tracking and stable
inversion: Enabling large feedforward tasks on a position-dependent
printer,” Mechatronics, vol. 38, pp. 76–90, 2016.
[11] T. Chen and B. A. Francis, Optimal Sampled-Data Control Systems.
London: Springer London, 1995.