印表機馬達系統之系統判別

國

立

交

通

大

學

電機與控制工程研究所

碩

士

論

文

印表機馬達系統之系統判別

System Identification of Motor System of Printer

研 究 生：洪朝瀛

指導教授：李福進 博士

印表機馬達系統之系統判別

System Identification of Motor System of Printer

研 究 生：洪朝瀛 Student：Chao-Ying Hung

指導教授：李福進 博士 Advisor：Dr. Fu-Ching Lee

國 立 交 通 大 學

電機與控制工程研究所

碩 士 論 文

A Thesis

Submitted to Department of Electrical and Control Engineering College of Electrical Engineering and Computer Science

National Chiao-Tung University in Partial Fulfillment of the Requirements

for the Degree of Master in

Electronics and Control Engineering

June 2005

Hsinchu, Taiwan, Republic of China

印表機馬達系統之系統判別

研究生：洪朝瀛

指導教授：李福進 博士

國立交通大學

電機與控制工程學系

碩士論文

摘 要

本篇論文是利用傳統的系統判別方法，最小平方估測法和最小

平方遞迴估測法加以改進，然後利用這些方法找出 HP inkjet 3325

彩色噴墨印表機直流噴墨馬達系統的最佳轉移函數，以便我們在設計

系統的控制器時，可以非常的準確。最後在本篇論文的這個系統判別

法則可以很容易的轉換成電腦輔助設計程式，並且可由實際系統的輸

入輸出資料來取得較佳的參數估測。

System Identification of Motor System of

Printer

Student: Chao-ying Hung Advisor: Dr. Fu-Ching Lee

Institute of Electrical and Control Engineering

National Chiao Tung University

Abstract

The purpose of this thesis is to utilize the traditional system identification

method, the least-square estimation method and the recursive least-square

estimation method to improve then utilize these methods to find out the

optimal transformation function of the DC motor system of HP 3325

color inkjet printer. An accurate model of the motor system can be used

to design a proper servo control system. At last, the estimation methods

in the thesis can be easily transformed into CAD programs, and we can

誌 謝

首先由衷的感謝我的指導老師 李福進博士，兩年來悉心的指導，使得本 篇論文能夠順利的完成，並且感謝在我的碩士生涯中陪著我度過的每個人，在此 致上我心中最崇高的敬意以及謝意。 另外也非常感謝一路走來，在背後支持我的父母和家人，在此也由衷感謝他 們，能夠讓我在沒有後顧之憂下，全心的完成學業，謝謝你們。

表 目 錄

表 3.3-1 誤差平方和的平均值比較表一 ... 29 表 3.3-2 誤差平方和的平均值比較表二 ... 30

圖 目 錄

圖 2.1-1 單一輸入和單一輸出的離散系統 ... 4 圖 2.1-2 多重輸入和多重輸出的離散系統 ... 6 圖 2.2-1 最小平方估測之幾何意義 ... 9 圖 2.3-1 實際輸出波形的局部放大圖 ... 16 圖 2.3-2 原始輸出波形與濾波後的輸出波形比較圖 ... 16 圖 3.1-1 模擬比較圖一 ... 19 圖 3.1-2 模擬比較圖二 ... 20 圖 3.1-3 模擬比較圖三 ... 20 圖 3.1-4 模擬比較圖四 ... 21 圖 3.1-5 模擬比較圖五 ... 21 圖 3.1-6 模擬比較圖六 ... 22 圖 3.1-7 二階系統的波德圖 ... 22 圖 3.1-8 模擬比較圖一 ... 24 圖 3.1-9 模擬比較圖二 ... 25 圖 3.1-10 模擬比較圖三 ... 25 圖 3.1-11 模擬比較圖四 ... 26 圖 3.1-12 模擬比較圖五 ... 26 圖 3.1-13 模擬比較圖六 ... 27 圖 3.1-14 三階系統的波德圖 ... 27

目 錄

中文摘要... i

英文摘要... ii

誌謝... iii

目錄... iv

表目錄... vi

圖目錄... vii

第一章 緒論

...

1

1.1 研究背景 ... 1 1.2 章節介紹 ... 1

第二章 定理及演算法

...

3

2.1 動態系統... 3 2.1.1 簡介 ... 3 2.1.2 線性差分方程式 ... 3 2.2 傳統系統判別方法... 6 2.2.1 最小平方估測法 ... 7 2.2.2 最小平方遞迴估測法 ... 10 2.3 估測法的改進... 15

第三章 實際例子模擬與比較

...

17

3.1 最小平方估測法的模擬... 17 3.2 最小平方遞迴估測法的模擬... 28

3.3 結果比較 ... 29

第四章 結論 ...

31

第一章 緒論

1.1 研究背景

要設計一個好的控制器，首先我們必須要找出系統的最佳轉移函數的數學模 式，在本論文中我們所要研究的就是如何以現有的傳統系統判別方法，最小平方 估測法和最小平方遞迴估測法加以改進，去找出系統最佳轉移函數的數學模式， 在本論文中我們是以 HP inkjet 3325 彩色噴墨印表機直流噴墨馬達系統的實際 例子來做研究，根據從工研院所量測到的輸入輸出資料來分析，再以傳統的系統

判別方法，最小平方估測法[1](Least-Squares Estimation Method)和最小平方

遞迴估測法[1](Recursive Least-Squares Estimation Method)和利用改進真實

系統的輸出值去改進的估測法，來加以計算以求得最佳的系統轉移函數的數學模 式，然後利用電腦輔助軟體(Matlab)去模擬並且和實際量測到的輸出資料做比 較，看是否吻合，來判別我們利用改進後的估測法，所求得的系統轉移函數的數 學模式，是否比未改進的估測法，所求得的系統轉移函數的數學模式，更接近真 實系統的轉移函數的數學模式，若我們所求得的系統轉移函數的數學模式越接近 真實的系統，那麼對我們要設計一個好的控制器將會有很大的幫助[2]。

1.2 章節介紹

第二章中我們將介紹最小平方估測法和最小平方遞迴估測法的理論基礎和 演算法的推導，以及利用改進真實系統的輸出值改進的估測法。第三章我們將會

利用第二章所講的方法去分析一個以 HP inkjet 3325 彩色噴墨印表機直流噴墨 馬達系統的實際例子（此平台是由工研院所架設的），所量測到的輸入輸出資料 去求得最佳的系統轉移函數的數學模式，並且利用電腦輔助軟體去模擬印表機馬 達系統之輸出波形，然後跟真實的系統所量測到的輸出波形互相作比較，以證明 改進後的估測法，所求得的參數估測值，比未改進前的估測法，所求得的參數估 測值，更接近真實的系統，第四章則為本論文之結論。

第二章 定理及演算法

在本章節裡，我們在 2.1 節中首先介紹動態系統及線性差分方程式，接著 2.2 節在介紹傳統系統判別方法，最小平方估測法及最小平方遞迴估測法，並舉 例說明，由於上述兩種方法都是利用最小平方誤差和來求得最佳解，加上印表機 馬達系統所量測的時間很短，我們所得到的資料有限，因此在 2.3 節中我們將介 紹如何利用改進真實系統的輸出值來改進估測法，使我們在使用相同數量的資料 下，卻能得到更接近真實系統的參數估測值。

2. 1 動態系統

2. 1. 1 簡介 動態系統可以被區分為兩種形式︰連續時間和離散時間。這兩種形式主要的 區別是在連續時間的系統中，訊號是連續的，而在離散時間的系統中，訊號也是 連續的，但是訊號是以離散的方式來表示。連續訊號的資料如果以適當的取樣比 率去取樣來儲存它的資料，那麼離散系統將會非常的接近連續系統[3]。 由於數位電腦的採用，所以系統判別技術絕大多數都是利用數位的方式來表 示，因此採用離散時間系統將會使我們更方便的去處理資料。在本論文中所提到 的最小平方估測法及最小平方遞迴估測法，也是建立在離散時間系統的基礎上。 2. 1. 2 線性差分方程式 首先，參考圖 2.1-1 的方塊圖，通常n th order− 差分方程式可表示成如下：

1 0 1 ( ) ( 1) _n ( ) ( ) ( 1) _n ( ) y k +a y k− + +… a y k−n =b u k +b u k− + +… b u k−n (2.1-1) or 1 0 ( ) ( ) ( n n j j j j y k a y k j b u k j = = +

∑

− =

∑

− (2.1-2) ) 其中y是系統的輸出， u 是系統的輸入， k 是整數時間 1, , n , 0 , 1, , a ⋅ ⋅ ⋅ a b b ⋅ ⋅ ⋅ b_n 是系統的參數。

Linear Discrete

System

u(k)

y(k)

圖 2.1-1：單一輸入和單一輸出（SISO）的離散系統 如果我們把 q 假設成位移運算子則 1 ( ) ( 1) q y k− = y k− 且 1 1 2 1 2 1 1 2 0 1 2 ( ) 1 ( ) n n n n A q a q a q a q B q b b q b q b q − − − − − − = + + + ⋅⋅⋅ + = + + + ⋅⋅⋅ + , − − ) (2.1-3) 則方程式(2.1-1)可被表示成下面的形式 1 1 ( ) ( ) ( ) ( A q− y k =B q− u k (2.1-4)

這是一個非常基本的數學模式，我們將會使用在系統判別（System Identification）上[4][5]。 仔細的觀察上面的式子，我們發現上面的推導過程可以很容易的去表示出差 分方程式和線性非時變的轉移函數之間的簡單關係。再來我們對方程式(2.1-1) 取 Z-轉換，並且假設初始狀況為零y k( )=u k( )=0，k< ，我們將會獲得下面的0 式子： 1 2 1 2 1 2 0 1 2 (1+ a z− + a z− + ⋅ ⋅ ⋅ + a z_n −n) ( )Y z = (b +b z− +b z− + ⋅ ⋅ ⋅ + b z_n −n)U z( ) 其中 z 是 Z-轉換的變數。 則轉移函數將可被定義為 1 2 0 1 2 1 2 1 2 ( ) ( ) ( ) 1 n n n n b b z b z b z Y z H z U z a z a z a z − − − − − − + + +⋅⋅⋅+ = = + + +⋅⋅⋅+ (2.1-5) or 1 1 ( ) ( ) ( ) B z H z A z − − = (2.1-6) 其中多項式 1 和 ( A z− ) B z( −1)定義在方程式(2.1-3)。 方程式（2.1-1）也可用來表示多重輸入和多重輸出的系統。考慮一個有 m 個輸入和 r 個輸出的的系統（如圖 2.1-2），然後我們將定義向量 和 如 下： ( ) U k Y k( ) 1( ) 1( ) ( ) ( ) ( ) ( ) m r u k y k U k Y k u k y k ⎡ ⎤ ⎡ ⎢ ⎥ ⎢ =_{⎢ ⎥} =_⎢ ⎢ ⎥ ⎢ ⎣ ⎦ ⎣ ⎤ ⎥ ⎥ ⎥⎦ 然後這個系統可以被表示成

1 0 ( ) ( ) ( ) n n j j Y k A Y k_j j B U k_j j = = +

∑

− =

∑

− (2.1-7) 其中A 和 j Bj分別是r r× 和r m× 的系統參數矩陣。 根據單一輸入和輸出的系統表示，我們把方程式(2.1-7)化簡成下列形式也 是可能的 1 1 ( ) ( ) ( ) ( A q− Y k =B q− U k) ) ) (2.1-8) 其中 1 和 ( A q− B q( −1 是以 1為變數的矩陣多項式，定義如下： q− 1 1 1 1 ( ) ( ) ( ) 1 ( ) _{0 1}( ) ( n n A q I A q A q n ) B q B B q B_n q − − − − =

+

+ +

= +

+ +

− − (2.1-9)

Linear Multivariable

Discrete System

u

₁

(k)

u

₂

(k)

u

_m

(k)

…

y

₁

(k)

y

₂

(k)

y

_r

(k)

…

圖 2.1-2：多重輸入和多重輸出(MIMO)的離散系統

2. 2 傳統系統判別方法

首先，考慮下面系統的數學模式 1 1 ( ) ( 1) _n ( ) ( 1) _m ( ) y k + a y k − + ⋅ ⋅ ⋅ + a y k − n = b u k − + ⋅ ⋅ ⋅ + b u k − m

其中 是系統的輸出， u 是系統的輸入，y a1,⋅ ⋅ ⋅, a n ,b1,⋅ ⋅ ⋅,bm 是系統的參 數。我可將其表示為 1 1 1 1 ( ) ( ) ( 1) ( ) ( 1) ( ) [ ( 1) ( ) ( 1) ( ) ] T n m n k m y k a y k a y k n b u k b u k m a a y k y k n u k u k m b b ϕ θ = − − − ⋅ ⋅ ⋅ − − + − + ⋅ ⋅ ⋅ + − ⎡ ⎤ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ = − − − − − − _{⎢ ⎥} ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎣ ⎦ 假設

[

]

T ( ) (1) (2) ( ) (1) ( ) ( ) T T Y k y y y k k k ϕ ϕ = ⎡ ⎤ ⎢ ⎥ Φ _{= ⎢ ⎥} ⎢ ⎥ ⎣ ⎦ (2.2-1) 寫成矩陣的形式為 Y = Φ (2.2-2) θ 其中Y是系統的輸出向量，Φ 是已知的系統資料向量，而θ則是未知且要估測的 參數向量。 2. 2. 1 最小平方估測法 最小平方估測 [9]是參數估測最常用的方法之一，它的基本想法是找出一組

參數使得實際所量測到的輸出與經由假設之數學模式所得到的輸出兩者的誤差 平方和為最小。因此，參數θ的選取是要使下面的函數最小化 2 1 1 ( , ) ( ( ) ( )) 2 k i V θ k y i y i = =

∑

− (2.2-3) 其中y i( )是實際量測所得到的輸出值，y i 是根據假設之數學模式所得到的輸出( ) 值。定義誤差值如下 ( )i y i( ) y i( ) y i( ) T( )i ε = − = −ϕ θ 則 2.2-3 式可寫成 2 2 1 1 1 1 ( , ) ( ) 2 2 2 k T i V θ k ε i E E E = =

∑

= = (2.2-4) 其中 E= − = − Φ Y Y Y θ 由於使 2.2-3 式函數最小化之參數θ滿足 ( , ) 0 T T V k Y θ θ θ _θ θ ₌ ∂ = −Φ + Φ Φ = ∂ T T Y θ ⇒ Φ Φ = Φ 若矩陣_{Φ Φ}T 是非奇異的矩陣，則 1 ( T ) Y θ = Φ Φ Φ (2.2-5) − T 2.2-5 式即為最小平方估測之演算法。

2 2

θ ϕ

1 1

θ ϕ

2

ϕ

1

ϕ

Y Y E 圖 2.2-1：最小平方估測之幾何意義 圖 2.2-1 以幾何的觀點說明最小平方估測的意義[9]；假設ϕ1 與ϕ2 為線性獨立之 向量，且 向量位於Y ϕ1_及ϕ2_{所組成的平面上，為使 E Y Y}= − 為最小，則 E必須 垂直於 所在的平面。因此參數Y θ的選取要使 剛好位於Y E的正交處；也就是Y 在此平面上的投影。 Example 考慮系統 1 2 ( ) ( 1) sin( ( )) ( 1) ( ) y t +ay t− =b u t +b u t− +e t 其中a= −0.8，b1 = −0.3， ， 為雜訊值，由於實際量測的資料通常含 有雜訊，所以在原來的系統上加上 項。 2 0.5 b = e t( ) ( ) e t 現在給予輸入訊號 u ，同時量測系統的輸出訊號y，給予之輸入訊號，必須 要有足夠的頻率成份，以確保系統各個模態均可被激發表現出來，此即所謂的持 續激發（persistent excitation）。在此我們假設 u 為一個頻率 100Hz 的方波

訊號。 Solution : 取 1000 組量測的資料(此 1000 組資料是利用 Matlab 軟體去計算出來的，e(t) 是一個隨機的雜訊，利用 Matlab 的“rand＂指令產生的)代入 2.2-5 式，我們可 得參數估測值為 1 2 0.8002 0.3043 0.4978 a b b θ ⎡ ⎤ _{⎡− ⎤} ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ = = −_{⎢ ⎥} ⎢ ⎥ ⎢_⎣ ⎥_⎦ ⎢ ⎥ ⎣ ⎦ 2. 2. 2 最小平方遞迴估測法 在這個章節中，我們將利用方程式(2.2-5)的最小平方估測的解來推導出最 小平方遞迴估測的演算法[1]。 假設向量方程(2.2-2)是由 組方程式所組成，我們以下面的方式來表示這 個方程式: m m m Y = Φ θ (2.2-6) 因此，我們將方程式(2.2-5)中的θ以θ( )m 取代 1 ( ) ( T ) m m m m m θ _{= Φ Φ} − ΦT_{Y (2.2-7)} 假設我們已經獲得一個新的方程式，則第m+ 方程式可表示如下 1 1 1 2 2 ( 1) ( 1) ( 1) _{n n}( 1) y m+ =θ ϕ m+ +θ ϕ m+ + +… θ ϕ m+ 定義

1 2 ( 1) [ ( 1), ( 1), , ( 1) T n m m m m ϕ + = ϕ + ϕ + ϕ + ] 則 (y m+ =1) ϕT(m+1)θ (2.2-8) 因此系統的第m+1方程式將可被寫成 1 m m Y ₊

=

Φ ₊₁θ (2.2-9) 其中 1 (1) ( ) ... ... ( 1) ( 1) m m y Y Y y m y m y m + ⎡ ⎤ ⎡ ⎤ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ = = ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ₊ ⎥ ₊ ⎣ ⎦ ⎣ ⎦ (2.2-10a) 1 1 1 1 (1) (1) ( ) ( ) ... ... ( 1) ( 1) ( 1) n m n m T n m m m m m ϕ ϕ ϕ ϕ ϕ ϕ ϕ + ⎡ ⎤ ⎡ ⎤ ⎢ ⎥ ⎢ _Φ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ Φ = = ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ _{+ +} ⎥ ⎢_{⎣ +} ⎥_⎦ ⎣ ⎦ (2.2-10b) 那麼新的最小平方估測器變為 1 1 1 1 1 (m 1) ( Tm m ) Tm m θ +

=

Φ Φ+ + − Φ +Y + (2.2-11) 從上式中(2.2-11)，我們很清楚的知道，如果要去求得新的最小平方估測器 (m 1 θ + ) ，我們必須去計算一個 n n× 矩陣的反矩陣，這樣反而會使我們的計算效 率變差。我們將去改進上面的缺失，利用舊的估測值θ( )m 去估測新的估測值 (m 1 θ + ) 而不需要去求反矩陣。在推導過程前，我們先介紹一下輔助定理: 輔助定理 : 假設A C, 和 A+BCD 是非奇異的方陣，然後 (A BCD+ )−1= A−1−A B C−1 ( −1+DA B−1 )−1DA−1

證明如下 證明:

[

]

1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 ( ) ( ) ( ) A BCD A A B C DA B DA I BCDA B C DA B DA BCDA B C DA B DA 1 I BCDA B I CDA B C DA B C DA B DA I BCDA BC C DA B C DA B DA I BCDA BCDA − − − − − − − − − − − − − − − − − − − − − − − − − − − − − − − ⎡ ⎤ + _⎣ − + _⎦ = + − + − + ⎡ ⎤ ⎡ ⎤ ⎡ ⎤ = + − _⎣ + _{⎦ ⎣} + _{⎦ ⎣} + _⎦ ⎡ ⎤ ⎡ ⎤ = + − _⎣ + _{⎦ ⎣} + _⎦ = + − 1 − I − = 故得證 首先，定義 1 ( ) ( Tm m) P m = Φ Φ − (2.2-12) 因此 1 1 1 ( 1) ( Tm m ) P m+ = Φ Φ+ + − 取代方程式(2.2-10)並且利用輔助定理，則我們可重寫P m( +1)如下： 1 1 1 ( 1) [ ( ) ( 1) ( 1)] ( ) ( ) ( 1)[1 ( 1) ( ) ( 1)] ( 1) ( T T T P m P m m m P m P m m m P m m m P m ϕ ϕ ϕ ϕ ϕ ϕ − − − + = + + + =

−

+ + + + ⋅ + ) (2.2-13) 從方程式(2.2-13)中我們可得知 1 1 ( 1) ( 1)[ ( 1) ( 1)] ( ) ( ) ( 1)[1 ( 1) ( ) ( 1)] ( 1) ( ) ( ) ( 1) ( 1) ( ) ( 1) [1 ( 1) ( ) ( 1)] ( 1) ( ) ( 1) ( 1 m m T m m T m m T T T T T m P m Y m y m P m Y P m m m P m m m P m Y P m m y m P m m m P m m m P m m y m θ ϕ ϕ ϕ ϕ ϕ ϕ ϕ ϕ ϕ ϕ − − + = + Φ + + + = Φ − + + + + ⋅ + Φ + + + − + ⋅ + + + ⋅ + + + ) ϕ 我們可以重新整理上式的最後兩項如下形式

1 1 ( ) ( 1)[1 ( 1) ( ) ( 1)] [1 ( 1) ( ) ( 1 ( 1) ( ) ( 1)] ( 1) ( ) ( 1)[1 ( 1) ( ) ( 1)] ( 1) T T T T P m m m P m m m P m m m P m m y m P m m m P m m y m ϕ ϕ ϕ ϕ ϕ ϕ ϕ ϕ ϕ ϕ − − + + + + ⋅ + + + − + + ⋅ + = + + + + ⋅ + ) 但是我們從方程式(2.2-7)和(2.2-12)得知 ( )_{m P m}( ) T_{m mY} θ = Φ 因此 (θ m+ ) 我們可以簡化成 1 1 ( 1) ( ) ( ) ( 1)[1 ( 1) ( ) ( 1)] [ ( 1) ( 1) ( )] T T m m P m m m P m m y m m m θ θ ϕ ϕ ϕ ϕ θ − + = + + + + + ⋅ + − + (2.2-14) 上面的結果顯示新的估測值將由舊的估測值加上一個校正值，而在校正項中的矩 陣 將可利用方程式(2.2-13)的遞迴公式來更新值，從(2.2-13)和(2.2-14) 這兩個公式中，我們可以很明顯的知道反矩陣的需求已經消除了，因此估測值 ( ) P m θ 更新的計算效率將會徹底改進。 接下來根據下面的遞迴演算法我們將會很容易的去實現最小平方遞迴估測 法： (m 1) ( )m (m 1) ( ) (P m m 1)[ (y m 1) T(m 1) ( ) θ + =θ +γ + ϕ + + −ϕ + θ m ] ) (2.2-15a) ( 1) ( ) ( 1) ( ) ( 1) T( 1) ( P m+

=

P m

−

γ m+ P m ϕ m+ ϕ m+ P m (2.2-15b) 其中 (γ m+ =1) 1 [1+ϕT(m+1) ( ) (P mϕ m+ ]1) 從上面式子我們可得知一個初始估測值 和相對應的 給定後且新的觀測 值可連續的獲得時，那麼我們將可連續不斷的去獲得新估測值 (0) θ P(0) θ。 最後，我們將對初始值θ(0)和P(0)的選擇來做討論，我們將在下面用兩個

實際的方法來討論： 1.取最先的 個資料且利用方程式(2.2-7)和(2.2-12)去求出 和 ，然 後利用所求的 和 去計算 m θ( )m P m( ) ( )m θ P m( ) θ(m+ 和1) P m( +1)。 2.假設θ(0)為任意數和P(0)=αI，其中α為很大的正數且I為單位矩陣，下面 將證明為何取此初始值來做計算。 證明： 利用方程式(2.2-15b)我們可得知 1 ( ) [ (0) T_{m m}] P m = P −

+

Φ Φ −1 (2.2-16) 因此 θ( )m =P m( )[Φ_mTY_m+P(0)−1θ(0)] (2.2-17) 為了使方程式(2.2-16)和方程式(2.2-12)相同，我們必須使P(0)−1趨近於零， 因此α趨近於無窮大，從方程式(2.2-17)中我們知道 可假設為任意數，因為 趨近於零。 (0) θ 1 (0) P − Example 題目同 2.2.1 範例 Solution : 若 (0) 100 , (0) 0P = I θ = ，我們可得參數估測值為： 1 2 0.7998 0.3033 0.4992 a b b θ ⎡ ⎤ _{⎡− ⎤} ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ = = −_{⎢ ⎥} ⎢ ⎥ ⎢_⎣ ⎥_⎦ ⎢ ⎥ ⎣ ⎦

2.3 估測法的改進

最小平方估測和最小平方遞迴估測，都是利用數學模式去計算誤差的平方和 的最小值，所得到的結果。在此我們將利用改進真實系統的輸出值 ，也就是 加上一個低通的濾波器，去濾掉比較高頻率的雜訊，使系統的相對誤差值 ( ) y i ( ) ( ) y i −y i 縮小且輸出的曲線更平滑，讓我們在計算參數估測值時能夠更精確， 使我們在設計控制器時，能夠有好的結果。 這個估測法的改進是改進真實系統的輸出值 ，首先，我們先取真實系統 的輸出值 每兩點之間的中點，然後連接每個中點，使輸出的波形能夠很接近 原始波形，這樣一來不但降低了誤差值，且可以使我們在計算時，得到更精確的 解，下面我們將推導改進後的輸出值變為多少？ ( ) y i ( ) y i 首先，假設 為系統的第i點輸出資料，然後其前後各點分別為 和 ，如圖 2.3-1 的實際輸出波形的局部放大圖所示，則我們可求得 ( ) y i y i( −1) ( 1) y i+ ( ) ( 1) 1 ( ) 2 ₂ ( ) ( 1) 1 ( ) 2 ₂ y i y i y i y i y i y i + − − = + + + = (2.3-1) 假設y i( )為改進後系統的第 點輸出資料，因為i y i( )和y i( −1)和 可看似三 角形的三頂點，且 ( 1) y i+ 1 ( ) 2 y i− 和 ₍ 1 2 y i+ ) 分別是兩邊的中點，所以y i( )落在 1 ( ) 2 y i− 和 ₍ 1 2 y i+ ) 連線的中點上，則由方程式(2.3-1)我們可得 ( ) ( 1) ( ) ( 1) 1 1 ( ) ( ) 2 2 2 2 ( ) 2 ( ) ( 1) ( 1) 2 4 y i y i y i y i y i y i y i y i y i y i 2 + − ₊ + − + + = = − + + = + + (2.3-2)

( ) y i ( 1) y i− ( 1) y i+ 圖 2.3-1:實際輸出波形的局部放大圖 為了證明我們改進後的輸出波形是否更平滑，接下來我們將上面方程式 (2.3-2)所求得的y i( )應用在實際的系統輸出資料中，則我們可從圖 2.3-2 中得 知輸出的波形如我們所想的變的比較平滑，且相對的誤差值也變的比較小。 圖 2.3-2:原始輸出波形與濾波後的輸出波形比較圖

第三章 實際例子模擬與比較

在本章節中，我們將會利用前一章所講到的演算法，來計算印表機之馬達系 統轉移函數的參數估測值，然後在去跟真實系統的輸出做比較，而在此我們所使 用的資料，是根據工研院所架設的印表機之馬達系統平台所量測到的。

3.1 最小平方估測法的模擬

在本節中，我們將利用最小平方估測演算法來模擬實際的例子(HP inkjet 3325 彩色噴墨印表機直流噴墨馬達系統)，然後我們將會在 3.3 節中去比較結 果，找出哪一個方法所求出的結果比較接近真實系統。因為我們所擷取的資料是 以離散的方式來儲存，所以在本節中，我們將會以離散時間系統來表示我們所要 求的系統，而且我們之前介紹的演算法都是以離散時間系統為基礎，故離散時間 系統的使用將會使我們能夠更方便的去計算。 首先，我們假設印表機之馬達系統是一個 2 階的離散時間系統[8]，我們定 義如下： 1 2 2 1 2 ( ) ( ) ( ) b z b Y z H z z a z a U z + = = + + 其中Y z( )是系統輸出，U z( )是系統輸入，a a b b1, 2, ,1 2是未知的參數估測值。 上式我們可以化簡成如下所示: 1 2 1 2 1 2 1 2 (1+a z− +a z− ) ( )Y z =(b z− +b z− ) ( )U z (3.1-1) 從 2.1 節中我們可得知方程式(3.1-1)將可轉換成差分方程式

1 2 1 2 ( ) ( 1) ( 2) ( 1) ( 2) y k +a y k− +a y k− =b u k− +b u k−

[

]

1 2 1 2 ( ) ( 1) ( 2) ( 1) ( 2) a a y k y k y k u k u k b b ϕ θ ⎡ ⎤ ⎢ ⎥ ⎢ ⎥ = − − − − − − ⎢ ⎥ ⎢ ⎥ ⎣ ⎦ (3.1-2) 再來我們將從馬達系統中所量測的輸入和輸出資料代入方程式(3.1-2)中去求出 和 ( ) y k ϕ 的值，然後將這些資料代入方程式(2.2-1)中，可求得 和 這兩 個矩陣，利用這兩個矩陣，我們代入最小平方估測演算法中(2.2-5)，我們可求 得 ( ) Y k Φ( )k 1 2 1 1 2 0.5342 0.3999 3.4524 3.5661 a a b b θ − ⎡ ⎤ ⎡ ⎤ ⎢ ⎥ ⎢₋ ⎥ ⎢ ⎥ ⎢ ⎥ = = ⎢ ⎥ ⎢− ⎥ ⎢ ⎥ ⎢ ⎥ ⎣ ⎦ ⎣ ⎦ 接著我們把 2-3 節所求的y i( )代入方程式(3.1-2)中去求出改進後的y k( )和ϕ 的 值，然後將這些值代入方程式(2.2-1)中，可求得Y k( )和Φ( )k 這兩個矩陣，利用 這兩個矩陣，我們代入最小平方估測演算法中(2.2-5)，我們可求得 1 2 2 1 2 0.5368 0.3972 3.4207 3.5406 a a b b θ − ⎡ ⎤ ⎡ ⎤ ⎢ ⎥ ⎢₋ ⎥ ⎢ ⎥ ⎢ ⎥ = = ⎢ ⎥ ⎢− ⎥ ⎢ ⎥ ⎢ ⎥ ⎣ ⎦ ⎣ ⎦ 在求得上面的參數估測值後，接下來我們將利用上面所求的參數估測值實際 的去模擬輸出波形及波德圖，然後檢測模擬出的波形是否符合印表機之馬達系統 的輸出波形，下面我們列出六個不同輸入狀況的模擬比較圖及它的波德圖

(3.1-7)，從圖 3.1-1 到圖 3.1-6 我們大致上可得知模擬出的波形都很吻合實際 系統的波形，只是我們無法從圖形中去得知θ₁和θ₂，哪一個參數估測值所模擬出 的波形比較接近真實系統，這個問題我們將在 3-3 節中利用模擬的輸出值相對實 際系統的輸出值的誤差平方和的平均值來比較哪一個比較接近真實系統。 而在下面的圖形中我們可看到original和 和 這幾個表示文 字，它們分別代表最小平方估測法和改進後的估測法和原始系統的輸出波形。 improved output 圖 3.1-1:模擬比較圖一

圖 3.1-2:模擬比較圖二

圖 3.1-4:模擬比較圖四

圖 3.1-6:模擬比較圖六

接著我們考慮一個 3 階的馬達系統之離散時間系統[8]，我們定義如下: 2 1 2 3 3 2 1 2 3 ( ) ( ) ( ) b z b z b Y z H z z a z a z a U z + + = = + + + 上式我們可以化簡成如下所示: 1 2 3 1 2 3 1 2 3 1 2 3 (1+a z− +a z− +a z− ) ( )Y z =(b z− +b z− +b z− ) ( )U z (3.1-3) 從 2.1 節中我們可得知方程式(3.1-3)將可轉換成差分方程式

[

]

1 2 3 1 2 3 1 2 3 1 2 3 ( ) ( 1) ( 2) ( 3) ( 1) ( 2) ( 3) ( ) ( 1) ( 2) ( 3) ( 1) ( 2) ( 3) y k a y k a y k a y k b u k b u k b u k a a a y k y k y k y k u k u k u k b b b ϕ θ + − + − + − = − + − + − ⎡ ⎤ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ = − − − − − − − − − _{⎢ ⎥} ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎣ ⎦ 再來如同上面 2 階的做法，我們可得到 1 2 3 1 1 2 3 0.5102 0.1761 0.2533 9.8145 16.7909 6.8623 a a a b b b θ − ⎡ ⎤ ⎡ ⎤ ⎢ ⎥ ⎢ ₋ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ − ⎥ =_{⎢ ⎥ ⎢ −}= ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ − ⎢ ⎥ ⎣ ⎦ ⎣ ⎦ ⎥ 1 2 3 2 1 2 3 0.5111 0.1774 0.2510 9.2984 15.7696 6.3635 a a a b b b θ − ⎡ ⎤ ⎡ ⎤ ⎢ ⎥ ⎢ ₋ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ − ⎥ =_{⎢ ⎥ ⎢}= _⎥ − ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ − ⎢ ⎥ ⎣ ⎦ ⎣ ⎦ 其中θ₁和θ₂分別代表最小平方估測值及改進後的最小平方估測值。 在求得上面的參數估測值後，接下來我們將利用上面所得到的參數估測值實 際的去模擬輸出波形及波德圖，然後檢測模擬出的波形是否符合印表機之馬達系 統的輸出波形，下面我們列出六個不同輸入狀況的模擬比較圖及它的波德圖 (3.1-14)，從圖 3.1-8 到圖 3.1-13 我們大致上可得知模擬出的波形都很吻合實

際系統的波形，只是我們無法得知參數估測值θ₁和θ₂哪一個比較接近真實系統， 這個問題我們將在 3-3 節中利用模擬的輸出值相對實際的輸出值的誤差平方和

的平均值來比較哪一個估測值的比較接近真實系統。

圖 3.1-9:模擬比較圖二

圖 3.1-11:模擬比較圖四

圖 3.1-13:模擬比較圖六

3.2 最小平方遞迴估測法的模擬

在本節中，我們將利用最小平方遞迴估測演算法來模擬實際的例子(HP inkjet 3325 彩色噴墨印表機直流噴墨馬達系統)，然後我們將會在 3.3 節中去 比較結果，找出哪一個方法所求出的結果比較接近真實系統。 如同 3.1 節的假設，我們將在此求出印表機之馬達系統的 2 階及 3 階的參數 估測值，首先，我們從 2 階的系統著手，將我們從馬達系統中所量測到的數據代 入最小平方遞迴估測演算法中，則我們可得到改進前的參數估測值θ₁及改進後的 參數估測值θ₂如下： 1 2 1 1 2 0.5321 0.4020 3.4500 3.5685 a a b b θ − ⎡ ⎤ ⎡ ⎤ ⎢ ⎥ ⎢₋ ⎥ ⎢ ⎥ ⎢ = = ⎢ ⎥ ⎢− ⎢ ⎥ ⎢ ⎥ ⎣ ⎦ ⎣ ⎦ ⎥ ⎥ 1 2 2 1 2 0.5340 0.4000 3.4407 3.5646 a a b b θ − ⎡ ⎤ ⎡ ⎤ ⎢ ⎥ ⎢₋ ⎥ ⎢ ⎥ ⎢ ⎥ = = ⎢ ⎥ ⎢− ⎥ ⎢ ⎥ ⎢ ⎥ ⎣ ⎦ ⎣ ⎦ 從θ₁和θ₂的數據上，我們可以知道最小平方遞迴估測法所求得的參數估測值 和最小平方估測法所求得的參數估測值誤差不大，所以我們可以很明顯的知道利 用最小平方遞迴估測法所求得的參數估測值，所模擬出的波形，也會吻合馬達系 統的輸出波形，故我們在此將不去做輸出波形的模擬比較。 接著我們考慮 3 階的系統，如同 2 階的做法，則我們可得到改進前的參數估 測值θ₁及改進後的參數估測值θ₂

1 2 3 1 1 2 3 0.5151 0.1739 0.2505 9.4984 15.9896 6.3849 a a a b b b θ − ⎡ ⎤ ⎡ ⎤ ⎢ ⎥ ⎢₋ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢− ⎥ =_{⎢ ⎥ ⎢−}= ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ − ⎢ ⎥ ⎣ ⎦ ⎣ ⎦ ⎥ 1 2 3 2 1 2 3 0.5162 0.1706 0.2528 9.6145 16.5009 6.7803 a a a b b b θ − ⎡ ⎤ ⎡ ⎤ ⎢ ⎥ ⎢₋ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢− ⎥ =_{⎢ ⎥ ⎢}= _⎥ − ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ − ⎢ ⎥ ⎣ ⎦ ⎣ ⎦ 3 階所求得的參數估測值和最小平方估測法求得的參數估測值接近，因此我們如 同 2 階一樣不去做輸出波形的模擬比較。 接著我們將在 3.3 節中討論未改進與改進後的最小平方遞迴估測法，哪一個 方法比較接近真實系統。

3.3 結果的比較

在本小節中，我們將去比較 3-1 節及 3-2 節所求出的參數估測值，去求出哪 一個比較接近真實系統，如下二表 3.3-1 及 3-3.2 所示： 表中 2 0 ( ( ) ( )) N i y i y i N ε = − =

∑

，其中y i( )是實際量測所得到的輸出值，y i 是根據( ) 假設之數學模式所得到的輸出值，N為系統的總資料數。 最小平方估測法： 表 3.3-1: 誤差平方和的平均值比較表一 系統階數 ε 未改進（ε₁） 改進後（ε ） ₂ 二階的系統 2.7582 1.8824 三階的系統 0.8889 0.8615

最小平方遞迴估測法： 表 3.3-1: 誤差平方和的平均值比較表二 系統階數 ε 未改進(ε₁) 改進後(ε₂) 二階的系統 1.6275 1.4293 三階的系統 1.1705 1.0516 從上面兩個比較表中，我們可以很明顯的知道，改進後的估測法，所求得的 參數估測值，所模擬出的輸出值相對於真實系統的輸出值的誤差平方和的平均值 比未改進的估測法，所求得的參數估測值，所模擬出的輸出值相對於真實系統的 輸出值的誤差平方和的平均值更小，因此我們知道改進後的估測法，所求得的參 數估測值比較接近真實系統。

第四章 結論

在印表機的馬達系統中，由於受限於噴墨頭的移動範圍不大，所以我們量測 資料的時間也相對的變短，量測資料時間的變短對我們所使用的系統判別方法將 會相當大的影響，因為我們所使用的最小平方估測法以及最小平方遞迴估測法， 都是需要使用越多的數據來計算則越精確。所以我們必須去找方法來改進，看能 不能在這樣的量測時間下，是否還能求出更接近真實系統的參數估測值，因此 2.3 節估測法的改進的使用是非常重要的，因為它降低了系統的誤差值，相對的 估測法的誤差平方和變小，因此我們從 3.1 節的模擬比較圖及 3.3 節的結果比較 表中，可以知道這個改進法讓我們得到更接近真實系統的參數估測值。 在本論文中，我們假設了 2 階及 3 階的兩個印表機之馬達系統，2 階系統的 假設，我們是根據自動控制系統書上[7]的馬達系統來假設的，而 3 階系統的假 設則是根據印表機的馬達系統的機械結構來分析的[6]，它是一個馬達系統加上 一個傳動系統，在這兩個系統中，我們必須去找出最接近真實系統的假設，則我 們從 3.1 節的模擬比較圖及 3.3 節的結果比較表中，我們可知道 3 階的馬達系統 比較接近真實的印表機之馬達系統。 最後，雖然我們可以求出較精確的印表機的馬達系統，但是在實際的系統中 有一個問題需要解決: 系統參數非固定: 例如:噴墨頭在向左與向右移動時，可能因與導桿間產生不同的摩擦力，以

及齒輪皮帶正反瞬間拉力不同，還有墨水在使用中的遞減，這些都會使系統參數

產生變化。

參 考 文 獻

[1] T. C. Hsia, University of California, Davis, System Identification: Least-Squares

Methods, Lexington Books, 1979.

[2] I. D. Landau, System Identification and Control Design, Prentice Hall, 1990. [3] Gene F. Franklin, et al, Digital Control of Dynamic Systems, Third Edition,

Addison-Wesley, 1998.

[4] Soderstrom Torsten, Stoica Petre, System Identification, Prentice Hall, 1988. [5] J. P. Norton, An Introduction to Identification, Academic Press, 1986. [6] R. Johansson, System Modeling and Identification, Prentice Hall, 1993.

[7] Benjamin C. KUO, Automatic Control Systems, Seventh Edition, Prentice Hall, 1995.

[8] M. Milanese and A. Vicino, “Optimal estimation theory for dynamic systems with set membership uncertainty : An overview ”, Automatica,Volume 27, Issue 6, pp. 997-1009, November 1991.

[9]張龍彬，「參數不確定系統判別法則與應用」，國立海洋大學電機工程學系 碩士論文，民國 88 年六月.