中 華 大 學

中 華 大 學 碩 士 論 文

題目：

^U

金屬濺鍍製程之批次控制

U

Run-by-Run Control of Metal Sputter Deposition

^U

系 所 別：機械與航太工程研究所 學號姓名：M09208007 郭子瑋 指導教授：陳俊宏 博士

中華民國 九十四 年 七 月

摘要

本論文中利用時間序列分析(Time Series Analysis)對金屬沈積率的歷史 資料建立濺鍍製程模型，然後將時間序列之製程模型轉換為狀態空間形 式，並利用延伸卡曼濾波器(Extended Kalman Filter)同時估計沈積率與更新 製程模型參數，據此改變輸入功率達到維持沈積率的目的。此外本文也提 供一個解決非等間距量測資料的預測靶材沈積率方法(Non-Fixed Interval Extended Kalman Filter)。文中最後將業界常使用的 Exponentially Weighted Moving Average (EWMA)與 Double Exponentially Weighted Moving Average (DEWMA) Controller 與本文提出的方法作比較，由比較結果可知， 無論是 沈積率的預測與控制，本文提出的方法皆有較好的效能表現。

關鍵詞：時間序列分析、延伸卡曼濾波器、濺鍍製程、沈積率

ABSTRACT

In this thesis, the time series model of the metal sputter deposition is constructed by history data and the Extended Kalman Filter is used to estimate the deposition rate and update the parameters of the time series model simultaneously. Then, the constant rate of deposition is accomplished by controlling the power based on the updated model. Furthermore, the thesis also proposes a method to estimate the deposition rate when the data is measured at non-fixed interval (Non-Fixed Interval Extended Kalman Filter). The simulation results demonstrate that the performance of the proposed method is higher than Exponentially Weighted Moving Average (EWMA) and Double Exponentially Weighted Moving Average (DEWMA) controllers that are adopted popularly in semiconductor processes.

Keyword: time series model, Extended Kalman Filter, deposition rate.

致謝

在這研究的兩年中，指導老師陳俊宏博士在論文上不辭辛勞的悉心指 導，並在程式上給予寶貴的意見，使本論文能順利完成，在此敬上最誠摯 的謝意。

在研究期間，感謝交通大學李安謙博士對論文之指導，以及陳文欽博 士、修平技術學院陳獻庚博士及南亞技術學院朱朝煌博士於口試時對於論 文的修改提供寶貴意見；感謝陳精一博士、許隆結博士及任貽明博士在研 究所期間課業及生活上的指點與照顧。

另外感謝實驗室的學長呂志男、學姐廖詩茵，在研究時期的幫忙，亦感 謝同學陳亭棋、張毓正、范牧樹與學弟李森正、劉佳原、廖乙安及葉祖銘 在課業與日常生活各方面的幫助。

最後，將本論文獻給我最摯愛的父母親及家人，謝謝你們不斷的給我支 持與鼓勵，使我順利完成學業。

目 錄

頁 次

中文摘要………Ⅰ 英文摘要………Ⅱ 誌謝………Ⅲ 目錄………Ⅳ 表目錄………Ⅵ 圖目錄………Ⅶ 符號說明………Ⅸ

第一章 緒論………1

1-1 研究動機………1

1-2 文獻回顧………1

1-3 研究方法與目的………5

第二章 直流電漿濺鍍介紹………7

第三章 Time Series 簡介………9

3-1 自我回歸過程(Autoregerssive Processes, AR)………12

3-2 移動平均過程(Moving Average Processes, MA)………12

3-3 自我回歸移動平均過程(Mixed Autoregerssive Moving Average Processes, ARIMA)………13

3-4 無定向型時間數列(Nonstationary Time Series)………13

3-5 季節性時間數列………14

3-6 Akaike’s Information Criterion (AIC)準則………16

3-7 時間序列與狀態空間方程式之轉換………16

第四章 Kalman Filter 簡介………18

4-1 狀態空間方程式簡介………18

4-2 Extended Kalman Filter………20

4-3 Non-fixed Interval Extended Kalman Filter………22

第五章 建立製程控制與補償機制………25

5-1 分析濺鍍製程資料………25

5-2 預測………30

5-3 製程模擬與控制………40

5-3-1 模擬方法………40

5-3-2 控制方法………41

第六章 結論………49

參考文獻………51

表目錄

表 1 Chamber 2 各靶材階次、係數與 AIC 值………28

表 2 Chamber 3 各靶材階次、係數與 AIC 值………29

表 3 DEWMA 與 Kalman Filter 之比較……… 30

表 4 DEWMA 權重………34

表 5 Chamber 2 兩種預測方法之 MSE 值………35

表 6 Chamber 3 兩種預測方法之 MSE 值………35

表 7 Age-Based DEWMA 之權重值………37

表 8 Chamber 2 兩種預測之 MSE 值………37

表 9 Chamber 3 兩種預測之 MSE 值………37

表 10 Chamber 2 等間距與非等間距預測結果比較………38

表 11 Chamber 3 等間距與非等間距預測結果比較………38

表 12 DEWMA Controller 模擬之權重………44

表 13 製程模擬預測沈積率之 MSE 值………44

表 14 製程模擬控制沈積率之 MSE 值………44

圖目錄

圖1 研究流程圖………6

圖2 直流電漿濺鍍機構造圖………8

圖3 模式建立流程………11

圖4 Kalman Filter運算流程圖………20

圖5 Non-Fixed Interval EKF 計算流程………24

圖6 Chamber 2 Target #1靶材沈積率歷史資料………26

圖7 擬合時間數列流程圖………27

圖8 兩種預測方法對Chamber 2 Target #1之沈積率預測結果比較圖……35

圖9 兩種預測方法對Chamber 3 Target #1之沈積率預測結果比較圖……36

圖10 兩種預測方法對Chamber 2 Target #1之沈積率預測結果比較圖……38

圖11 兩種預測方法對Chamber 3 Target #1之沈積率預測結果比較圖……39

圖12 模擬產生沈積率與實際沈積率比較圖………41

圖13 濺鍍控制架構圖………42

圖14 第一個靶材之模擬沈積率預測結果………45

圖15 第一個靶材之模擬沈積率預測結果，細部圖………45

圖16 第一個靶材之模擬沈積率控制結果………46

圖17 第一個靶材之模擬沈積率控制結果，細部圖………46

圖18 第二個靶材之模擬沈積率預測結果………47

圖19 第二個靶材之模擬沈積率預測結果，細部圖………47

圖20 第二個靶材之模擬沈積率控制結果………48

圖21 第二個靶材之模擬沈積率控制結果，細部圖………48

符號說明

A :狀態轉移參數 ( )

a k ：在 k 時期之白噪音 B ：後移算子

b ：斜率 C ：常數

d ：未量測之批次數 ( )

d k ：在 k 時期與 k-1 時期未量測之批次數 E ：期望值運算子

e ：輸出誤差 e ：預測誤差 p

F ：線性化狀態轉移矩陣

f ：非線性狀態函數 H ：線性化輸出矩陣

h ：非線性輸出函數 I ：單位矩陣

J ：目標函數

j ：大於零之正整數

( )k

K ：卡曼增益向量

k ：離散時間點 m ：大於零之正整數

( )k

P ：狀態共變異矩陣

ˆ ( )

p k

：漂移項估計值

( ) k

Q ：系統干擾之共變異矩陣

( ) k

R ：量測干擾之共變異矩陣

T ：沈積率目標值

( )

u k

：控制輸入

( )

v k

：量測干擾

( )

W k

：差分數列 ( )

W k ：減去平均值後之差分數列

1

,

2

w w

：權重

( )

w k

：製程干擾

( ) k

X ：系統狀態向量

( )

Y k

：沈積率觀測值數列 ( )

Y k ：沈積率觀測值減去常數後之數列

α

：截距

β

：轉移係數

( , ) k j

δ

：步階函數 Γ ：干擾輸入向量

( ) k

ε

：製程模型干擾項 Φ ：自我回歸係數矩陣

Φ

P：季節性自我回歸過程 P 階之係數

φ

p：自我回歸過程 p 階之係數

Θ ：季節性移動平均過程 Q 階之係數 Q

θ

q： 移動平均過程 q 階之係數

∇ ：差分運算子

上標：

D：季節性差分次數

d：差分次數

( )m ：排列順序

s：季節性數列之週期 下標：

P：季節性自我回歸過程之階次

p：自我回歸過程之階次

Q：季節性移動平均過程之階次

q：移動平均過程之階次

第一章 緒論

1-1 研究動機

隨著國內製造業的逐漸升級，電子化整合製造是工廠趨向高效能 操作管理的重要發展趨勢，當工廠進入量產以後，提昇製造良率對於 量產為主之製造系統^T如IC、LCD、或光碟片等 3C產業^T至為重要。但目 前濺鍍製程在上述產業中因設備與技術不斷提升，對於產品的品質要 求也相對的嚴苛，因此各產業無不尋求更佳的製程整合與控制系統，

來因應時代潮流趨勢。

為了解決上述問題，並進一步增進產品品質，透過先進製程監控 早期的預知診斷可能發生的製程異常偏離，並進一步的調整製程機台 參數，使得製程機台回歸正常狀態，藉由此一機制不僅可以將製程控 制在更穩定的狀態下，更可以正確的掌握機台停機維修時間、減少因 無預警停機而造成的在製品損失、有效提高人力運用效率、減少產品 品質變異、提升製程與設備有效操作率。

1-2 文獻回顧

製程控制方面Kalman[1]由狀態方程推導出離散遞回式Kalman Filter，此後便廣泛應用在各領域之狀態與參數估計上。Sachs et al.[2]

對於超大型積體電路(VLSI)製造廠提出一個包含三個模組的製程控制 系統，其模組分別為Flexible Recipe Generator、Run by Run Controller 和Real Time Controller，並以LPCVD製程驗證此製程控制系統之效 能。Spanos et al.[3]應用時間序列(Time Series)分析電漿蝕刻(Plasma Etching)製程參數，再利用T^P2P統計分析(Hotelling’s T^P2P Statistic)對先前建 立的時間序列模型所求得之殘差(a^Bt^B, White Noise)來做分析。經上述步 驟處理之資料可由T^P2P Control Chart來判讀製程發生的變異參數為何，以 及做預防保養(Predicted Maintain, PM)的決策時機。Butler與Stefani[4]

使用反應曲面法(Response Surface Model, RSM)找出輸入與輸出之數 學模型，針對量測值與數學模型輸出值之誤差，利用(Predictor-Collector Controller, PCC)的控制方法對閘極蝕刻製程產出良率做改善，經由實 驗證明使用PCC控制器後，能夠減少 36%的平均蝕刻率(Mean Etch Rate) 之標準誤差。Sachs et al.[5]提出一個完整的Run by Run架構。其中在控 制部分，先由Control Chart判別變異之類型後，啟用不同程度之修正模 式。修正模式分為量測值之短期大幅度的跳動(Large Shift)以及長期製 程漂移(Drift)現象，而處理方法分別為貝氏近似法(Bayesian Approach) 及EWMA。Boning et al.[6]討論在CMP製程中四個輸入變數(Speed, Pressure, Force, Profile)及兩個輸出變數(RemovalRate, Non-Uniformity) 的MIMO(Multiple Input-Multiple Output)的問題，其利用向量式線性回

歸及EWMA Controller補償誤差值方式來提升良率。Krauss et al. [7]對 於聚合物沈積製程，利用多個Kalman Filter組合之Interactive Multiple Model(IMM)估計製程參數；經模擬結果可知，當製程參數隨製程變動 時，IMM方法可以減少輸出的變異。Smith et al.[8]探討在CMP製程中 研磨率與不均勻度會產生漂移現象，利用類神經網路計算當時期漂移 與雜訊值，輸出修正的EWMA權值，再利用EWMA控制器修正漂移現 象。Sattler et al.[9]利用經驗公式建立金屬濺鍍製程模型，並使用EWMA 控制器來達到降低成本(33%)，提高產能(3.6%)之目的。Smith et al.[10]

討論在濺鍍過程中，沈積率(Deposition Rate)會隨著靶材使用的時間而 降 低 。 故 應 用 能 夠 有 效 消 除 製 程 漂 移 (drift) 的 控 制 方 法 (Predicted Corrector Controller, PCC)，來調整沈積時間，使沈積厚度維持在目標 值之間。Krauss 與Kamen [11] 在大面積平版之薄膜塗佈製程上，利用 Kalman Filter估計薄膜之黏滯係數，並且配合Run by Run製程控制，設 定是否更新製程模型係數之黏滯係數門檻值，以維持薄膜的厚度。Guo et al.[12]提出一種能夠偵測並使用修正權值使發生大位移製程能快速 回到穩定狀態的控制器(Self-Tuning EWMA Controller)。經由蒙地卡羅 模擬大位移及漂移狀況驗證後，證明能有效處理製程中所發生的位移 及漂移現象。Guo et al.[13] 對於控制確定性漂移，提出一套權重函數 法來求得控制器之最佳權值，其原理為降低暫態控制值，並且不增加

製程變異條件下，計算其目標函數製成圖表。再經由決定控制與變異 兩者比例，即可決定控制器最佳參數。Chen與Kuo[14] 在CMP製程中 考 慮 批 次 量 測 時 有 中 間 間 隔 一 段 時 間 才 批 次 量 測 狀 況 ， 在 原 有 DEWMA修正模型中再加入量測時間之變化量來修正移除率(Removal Rate)的變化率(trend)，使模型在修正漂移現象時更加準確。Tseng et al.[15]提及在過去使用EWMA模式時，權重的值由經驗來決定，而且 通常是定值，其次對於修正製程Shift現象需要較長的時間。所以針對 此問題，在原本固定權值外，加入隨時間指數遞減之權值，能夠在 5~15 run 內 修 正 回 目 標 值 。 Chemali et al.[16] 運 用 PROLITH 模 擬 曝 光 機 (Stepper)製程，利用Kalman Filter估計製程中參數，控制焦距(Focus)與 能量(Dose) 兩個變數，使CD(Critical Dimension)與SWA(Sidewall Angle) 控制在目標值之間。Corigliano與Mariani [17]討論Extend Kalman Filter 在狀態與參數估計時，對於單自由度與多自由度振動系統中彈簧軟化 (Softening)情況時，會出現系統狀態估計是準確的，而參數的估計與更 新模型卻造成系統模型與實際系統的偏差現象。

1-3 研究目的與方法

在於過去濺鍍製程大多是單機操作，與上下游製程並無關連，主要 目的在使濺鍍製程能達到^TRun By Run^T製程控制，並且將來能透過整廠 的先進製程控制^T(Advance Process Control, APC)^T來整合上下游的製程。

所進行的步驟為控制製程之輸入參數，減少濺鍍過程中製程變異程 度，達到產能及良率提昇。在濺鍍過程中靶材的消耗，反應腔體內原 子堆積或製程中內外部干擾的影響下，都會造成產品不良率上升。故 將以及時間序列分析^T(Time Series Analysis)^T來分析輸入參數、量測數據 及干擾之間的關係，再使用卡曼濾波器^T(Kalman Filter)^T估計製程中靶材 沈積率以及製程參數，並即時更新製程模型，最後加以穩定控制製程 使沈積率維持在目標值之間。

研究流程如圖 1 所示，起始步驟為收集相關資料及金屬濺鍍製程的 瞭解。第二步為分析濺鍍製程資料，第三步驟為針對金屬濺鍍性質建 立控制及補償機制，第四步驟為利用第三步驟所建立之控制法則進行 預測沈積率以及模擬製程控制，最後為結果與討論。

圖1 研究流程圖 收集資料及瞭解金屬濺鍍製程

分析濺鍍製程資料

建立製程控制與補償機制

模擬驗證

結果與討論

第二章直流電漿濺鍍介紹

直流電漿濺鍍機(DC Plasma Sputter)已經普遍的使用在半導體製造業，

例如鋁合金金屬導線沈積。這項製程技術也廣泛的應用在光碟片與LCD產 業中，直流電漿濺鍍是在真空艙體中，將兩極通上直流電，靶材置於陰極，

基材置於陽極，電漿產生於兩極間，兩電極中間部分氣體離子化後加速撞 擊置於陰極的靶材，使靶材中的原子脫離靶材，以自由碰撞或擴散的形式 沈積在基材上，如圖2所示。但是本製程只限於用在靶材為導電體上，因此 多用在金屬的濺鍍上。而為了提高電漿的濺鍍速率，一般都在靶材的背面 加上磁場。藉由磁場及電場的交互作用可使二次電子集中在靶材表面附 近，呈螺旋狀運動，使此一區域氣體分子的電離率提高，增加電漿密度，

可以擊出更多的靶材原子。整個沈積過程大概可以分為以下步驟[18]:

1. 電漿內氣體被離子化，將脫離電漿並往負極版移動 2. 經加速之離子撞擊(Bombard)靶材，將靶材原子撞出 3. 被撞擊出的原子經過電漿，傳遞到基材表面上 4. 原子吸附在基材上，進行沈積

圖2 直流電漿濺鍍機構造圖[18]

由以上的製程步驟可知濺鍍製程中濺鍍氣體的壓力、沈積時間、製 程溫度與輸入功率的大小與薄膜合金的穩定度，薄膜的厚度與平坦度有 著相當大的關係。

基材 基材

基座

第三章Time Series 簡介

在本論文中，將利用時間序列分析建立濺鍍製程中沈積率(Deposition Rate)與功率(Power)之關係，以及製程中干擾(Noise)的型態。之後將時間序 列數學模型轉換為狀態空間方程式；使用Kalman Filter 估計在每一個批次 之沈積率與係數狀態，並且預測下一個批次之沈積率。以下小節中將會介 紹時間序列與分析製程資料步驟。

所謂時間序列意指隨時間連續觀察所產生有順序的觀測值之集合。或 對於製程而言，時間序列即為隨時間連續觀測之量測值之集合。在處理時 間序列資料時，傳統計量模式所使用的迴歸方法，雖然簡便，但對序列資 料之自我相關以及誤差項的處理，總不夠嚴謹。而Box 與 Jenkins 於 1970 年代初期所提出的 ARIMA 模式，即可針對前述缺失作一改善[19-21]。

ARIMA 模式為利用變數本身之落後值，以及當期與過去隨機誤差項，加 權後以解釋變數本身之變動。Box 與 Jenkins 所提出之模式建構程序為一 種試誤遞迴過程(Trial and Error Iterative Process)，如圖 3 所示，以下將針 對此建構模式程序加以說明。

1. 資料蒐集

針對需進行時間序列分析之資料，蒐集整理成歷史性時間序列 資料。或在製程方面，現場製程機台與量測機台的資料，蒐集整理

成歷史性時間序列資料。

2. 製程模式的建立

在進行時間序列分析之前，必須確定時間序列是否為恆定 (Stationary)，如此才可以以一個固定係數的方程式來進行估計及預 測。ARIMA(p,d,q)模式的建構第一步須先決定差分(Differences)的 階次d ；接著再決定自我迴歸之落差期數，若落差p 期，則記為 AR(p)；同時也需決定移動平均期數，若落差q 期，則記為MA(q)。

將蒐集整理完成之時間序列資料予以繪圖，判定時間序列資料是 否恆定，若否，則加以差分直至數列恆定，並據以決定模式階次。

此外，再計算序列資料之自我相關函數(Autocorrelation Function, ACF)、偏自我相關函數(Partial Autocorrelation Function, PACF)，以 及延伸之自我相關函數(Extend Autocorrelation Function, EACF)，來 決定數列之落差期數。

若時間序列資料出現週期性時，則需對資料做進一步的處理。

亦即在原本模式中再加入週期性影響因子，此時模式可稱為季節 性(或週期性)時間序列(Season ARIMA, SARIMA)。同樣的，經過 上述鑑定模式之步驟後，可獲得更精準的數學模型。

3. 製程模式之參數估計

時間序列模型其變數必須不會隨時間的經過而發散、呈現不穩

定的特性，如此才能用固定係數的方式加以預測或估計。換句話 說，任何屬時間序列的資料必須具有恆定的特性，模型的分析與估 計才具有意義。

對於模式裡的參數

θ

、

φ

可利用最大概似函數(Max Likelihood Function)來求得。若時間序列的模式為非線性時，不能利用一般 的最大概似函數來估計模式中的參數，必須採用非線性的反覆求 解過程來估計模式中的參數。

4. 模式檢核

完成模式參數估計後，必須對模式進行適切性檢核。診斷時，

乃利用原始實際值與模式預測值之間的差距，即殘差值(Residual Value)，進行序列殘差值的自我相關程度檢定。

圖3 模式建立流程

蒐集資料

暫定模式

估計參數

診斷模式是否正 確 否

是

狀態空間模型

3-1 自我回歸過程(Autoregerssive Processes, AR)

AR 過程為當期觀測值與所有前期觀測值之回歸組成，即　 　 　 　 Y k( )= +C

φ

₁Y k( − +1)

φ

₂Y k( − +2) "+

φ

_pY k( − p)+a k( ) (3.1) 式中Y(k)表示觀測值，C 為常數，

φ

為權值，a(k)為白噪音(White noise)，

下標表示時序。故(3.1)式可稱為 p 階自我回歸過程(Autoregerssive Process of Order p, AR(p))。藉由後移運算子 B(Backward Operator)來表 示(3.1)式可得

1 2

1 2

( ) ( _p ^p) ( ) ( )

Y k = +C

φ

B +

φ

B +"+

φ

B Y k +a k (3.2) 令

1 2

1 2

( ) 1 ^p

p B B B pB

φ

= −

φ

−

φ

− −"

φ

(3.3)

則(3.2)式可表示為

( ) ( ) ( )

p B Y k C a k

φ

= + (3.4) 故(3.4)式即為 AR(p)之通式。

3-2 移動平均過程(Moving Average Processes, MA)

MA 過程其意義為當期觀測值是由當期及過去干擾值之線性組合，

即

1 2

( ) ( ) ( 1) ( 2) ( ) ( ) ( )

q

q

Y k a k a k a k a k q

B a k

µ θ θ θ

µ θ

= + − − − − − −

= +

"

(3.5)

模式中之係數 θ 亦稱為震動影響或記憶函數(Shock-Effect or Memory Function)，此即表示干擾影響會持續 q+1 個時期後消失。(3.5)式稱之為 q 階移動平均過程(Moving Average Process of Order q, MA(q))。(3.5)式使 用移動平均這個名稱僅說明 MA(q)過程為數個 a(k)之移動線性組合而 已，並非真正的移動平均，此乃因權數之和不等於1。

3-3 自我回歸移動平均過程(Mixed Autoregressive Moving Average Process, ARMA)

將上述 AR、MA 合併後可得一

(p,q)

階混合自我回歸移動平均過程 (Mixed Autoregressive Moving Average Process of Order(p,q))，其形式為

( ) ( ) ( ) ( )

p B Y k C q B a k

φ

= +

θ

(3.6)

事實上前述之AR(p)與 MA(q)皆為 ARMA 之特例，在擬合資料時 ARMA 模式會比單獨使用 AR 或 MA 模式來的精簡。並且 ARMA 包含了 AR 之平穩型與 MA 之可逆轉(Invertibility)兩者之特性，即當

φ

_p( ) 0B = 之根 落於單位圓外時，ARMA(p,q)為平穩型。當

θ

_q( ) 0B = 之根落於單位圓外 時，ARMA(p,q)則具有可逆轉性。

3-4 無定向型時間序列(Nonstationary Time Series)

前述之數列皆在恆定狀態下所求得，所謂之恆定即為時間序列中有

恆定的狀態，此時稱此數列為無定向數列。也就是說數列之平衡水準會 隨時間變動。為了能將數列做分析，必須對無定向數列做差分後，期使 數列能從Non-stationary 轉變為 Stationary 便於分析。在此定義差分運算 子(Difference Operator)為

( ) ( ) ( 1) (1 ) ( ) Y k Y k Y k B Y k

∇ = − − = − (3.7)

故經差分後之時間序列稱之為(p,d,q)階之整合自我回歸移動平均模式 (Autoregressive Integrated Moving Average Model of Order(p,d,q))，或簡稱 為ARIMA(p,d,q)，其模式之形式可表示為

( ) ^d( ( ) ) ( ) ( ) or ( ) ( ) ( ) ( )

p B Y k q B a k p B W k q B a k

φ

∇ −

µ

=

θ φ

=

θ

(3.8)

若在(3.8)式中增加常數項 C，則可變化成更一般化之整合自我回歸移動 平均過程，即

( ) ( ) ( ) ( )

p B W k C q B a k

φ

= +

θ

(3.9)

若在模式(3.9)式中令 C = 0 則此模式所代表之數列將具有隨機性趨勢 (Stochastic Trends)，即數列會在不同水準與不同斜率上隨機變化。若 C ≠ 0，則代表數列具有一個固定性 d 階多項式趨勢。而此時差分數列之平 均值為

1 2

( ( )) ( ^d ( )) _w /(1 _p)

E W k = E ∇ Y k =

µ

=C − − − −

φ φ

"

φ

(3.10)

3-5 季節性時間序列

許多的製程資料大多含有季節性(週期性)因素。例如：每道製程電 流變化量、腔體壓力、靶材的壽命…等。這些資料都有一共同特性，相 差固定時間間隔的觀測值，彼此間有高度的相關性。滿足上述特性之數 列可稱之為季節性時間序列。對原始數列做季節性(週期性)修正，其目 的在於消除季節落差(Seasonal Lag)的現象。季節性 ARMA 模式或 SARMA(P,Q)s 之形式為

( ) ( )^s ( ) ( )^s

P B Y k Q B a k

Φ = Θ (3.11)

其中Φ_P( )B^s 與Θ_Q( )B^s 分別表示為P 與 Q 階之B 的多項式，若為無定向^s 數列時則需對數列差分，則模式改寫為

( )(1^{s s D}) ( ) ( ) ( )^s

P B B Y k Q B a k

Φ − = Θ (3.12) (3.12)式中(1−B^{s D}) 為季節性差分，模式可記為ARIMA(P,D,Q)s，s 為季 節性資料之週期，例如以月為單位s=12、以季為單位 s=4。在實務應用 上，除了相鄰固定時間間隔之觀測值有相關外，連續的數個觀測值間也 可能有很高的自相關性。因此可用更一般化的模式，來描述此種數列的 特性。此模式稱為(p,d,q)*(P,D,Q)s 季節性相乘 ARIMA 模式，記為 SARIMA(p,d,q)(P,D,Q)s

( ) ( )(1^s ) (1^{d s D}) ( ) ( ) ( ) ( )^s

p B P B B B Y k q B Q B a k

φ

Φ − − =

θ

Θ (3.13)

其中

2

1 2

2

1 2

2

1 2

2

1 2

( ) 1 ( ) 1 ( ) 1

( ) 1

p

p p

s s s Ps

P P

q

q q

s s s Qs

Q Q

B B B B

B B B B

B B B B

B B B B

φ φ φ φ

θ θ θ θ

= − − − −

Φ = − Φ − Φ − − Φ

= − − − −

Θ = − Θ − Θ − − Θ

"

"

"

"

(3.14)

3-6 Akaike’s Information Criterion (AIC)準則

假設一組資料可以一個含有M 個參數的統計模式來擬合，作為評估 模式擬合的品質，Akaike 於 1973 年提出一種判定準則，其定義為：

2ln( ) 2

AIC = − ML + M (3.15) 對於ARMA 模式且含有 n 個有效觀測值之數列，其對數形式之概似函 數為

2

2

ln ln 2 1 ( , , ) 2 ^a 2 _a

L n

πσ

S

φ θ

c

= − −

σ

(3.16) 對參數

φ θ

, ,c 與

σ

_a²取(3.16)式為最大，則可得

ˆ2

ln ln (1 ln 2 ) 2 ^a 2

n n

L= −

σ

− +

π

(3.17)

上式中的第二項為常數，故AIC 可簡化為

ˆ2

ln _a 2

AIC =n

σ

+ M (3.18)

3-7 時間序列與狀態空間方程式之轉換

在此小節將說明如何將時間序列模型轉換為狀態空間方程式。現假

設時間序列模型為 ARIMA(p,d,q)過程，若令m=max(p+d q, + ，則1) ARIMA(p,d,q)模型可寫為[22]

( ) ( )B Y k ( ) ( )B a k

ϕ

=

θ

(3.19) 或

1 1

( ) ( 1)

_m

( ) ( ) ( 1)

_m

( 1 )

Y k = ϕ Y k − + " + ϕ Y k − m + a k − θ a k − − − " θ a k − + m

(3.20) 令狀態X 表示為

( )

( ) ( 1)

1

( )

m

m m

X k = ϕ Y k − − θ

₋

a k

(3.21)

( ) ( 1)

( ) ( 1) ( 1) 1 ( )

j j

m j

X k =

ϕ

Y k− + X ⁺ k − −

θ

₋a k (3.22) 其中 j m< ，

θ

₀

= − 1

，則狀態向量X

( ) ( k = X

⁽¹⁾

( ), , k " X

^{( )m}

( )) k

^T可表示為 下面的形式

1

1

1

1

1 0 1

( ) ( 1) ( )

0 1

0 0

m

m m

k k a k

ϕ

θ ϕ

ϕ θ

−

−

⎡ ⎤ ⎡ ⎤

⎢ ⎥ ⎢ − ⎥

⎢ ⎥ ⎢ ⎥

= − +

⎢ ⎥ ⎢ ⎥

⎢ ⎥ ⎢ − ⎥

⎣ ⎦ ⎣ ⎦

X X

"

# # % #

" #

"

(3.23)

而輸出方程可表示為

[ ]

( ) 1 0 0 ( )

Y k = " X k

(3.24) (3.23)式與(3.24)式即為時間序列模型(3.20)式轉換成狀態空間方程式之 形式，關於狀態空間方程式之介紹，會在下一章中說明。

第四章Kalman Filter 簡介

Kalman Filter 最早是 1960 年 Kalman [1]所提出，其利用狀態空間方程式 (State Space Equation)之特性所發展出的遞迴估計系統狀態方法。整體而言 卡曼濾波器具有以下之特點[23]：

1.卡曼濾波器之解是一個適合在電腦上計算的遞回方程式。

2.對於數據可以逐一的即時處理，並基於目前與前一刻計算出的狀 態估計值運用遞迴方式算出下一個狀態的估計值。

3.系統方程可以是時變的動態系統，輸入之訊號或數據之型態可以 是非平穩的。

4.可應用於各種線性(Linear)或非線性(Non-Linear)問題。

以下將介紹狀態空間方程式及 Extended Kalman Filter(EKF)之架構。

4-1 狀態空間方程式簡介

一個線性離散時間系統之狀態空間方程式可由狀態方程式與輸出 方程式之組合來表示，即[23]

(k + =1) (k +1, ) ( )k k + (k +1, ) ( )k k

X F X Γ w (4.1) (k + =1) (k +1) (k + +1) (k+1)

Y H X v (4.2) (4.1)式為狀態方程式，代表系統中狀態的變化，(4.2)式為輸出方程式，

代表由系統狀態產生變化後系統的輸出。其中 X 代表系統狀態向量

(n× ，Y 為系統輸出向量(1) m× ，F 為 k 到 k+1 時期之狀態轉移矩陣1) (n n× ，Γ 為干擾輸入矩陣() n n× ，H 為觀測矩陣() m n× ；w 與 v 分別) 為系統干擾向量 (n× 以及量測干擾向量(1) m× ，並且為 White-Noise，1) 且與系統狀態X 三者互相獨立，並且有以下性質：

[ ( ) ^T( )] ( , ) 0, , 0

E w k w j =Q

δ

k j Q≥ ∀k j≥ (4.3) [ ( ) ^T( )] ( , ), 0, , 0

E v k v j =R

δ

k j R≥ ∀k j≥ (4.4) [ ( ) ^T( 1)] 0, , 0

E w k v j+ = ∀k j≥ (4.5) [ (0) ^T( )] 0, [ (0) ^T( 1)] 0 0

E X w k = E X v k + = ∀ ≥k (4.6) 則與(4.1)式與(4.2)式相對應之卡曼濾波器形式為

ˆ(k + =1) ˆ(k +1 )k + (k +1)[ (k + −1) ˆ(k +1 )]k

X X K Y Y (4.7) ˆ(k +1 )k = (k +1, ) ( )k ˆ k

X F X (4.8) ˆ(k +1 )k = (k +1) (ˆ k +1 )k

Y H X (4.9)

K ( k + = 1) P ( k + 1 ) k H

^T

( k + 1)[ ( H k + 1) ( P k + 1 ) k H

^T

( k + + 1) R ( k + 1)]

⁻¹

(4.10)

( k + 1 ) k = ( k + 1, ) ( ) k k

^T

( k + 1, ) k + ( k + 1, ) ( ) k k

^T

( k + 1, ) k

P F P F Γ Q Γ

(4.11)

(k + = −1) [ (k+1) (k +1)] (k +1 )k

P I K H P (4.12) 初始值設定為

ˆ (0) = E [ (0)] = (0), (0) = var (0) ∀ ≥ k 0

X X X P X

(4.13) 式中K 為卡曼增益(Kalman Gain)，遞迴計算流程如下圖 4，利用已給定

之初始值預測下一個時點狀態變異矩陣P(k+1|k)，接著計算卡曼增益，

利用這些已知數據去預測下一個時點(k+1)之系統狀態 ˆ (X k +1 )k 及輸出 觀測值 ˆ (Y k+1 )k 。此時新的觀測值 Y(k+1)進入後，就可以更新目前在 k+1 時點的狀態估計向量

X ˆ ( k + 1)

以及狀態變異矩陣 (P k+1)。如此可計 算下一個狀態變異矩陣P(k+2|k+1)之預測值，依次類推。

圖 4 Kalman Filter 運算流程圖[23]

4-2 Extended Kalman Filter

Extended Kalman Filter 是將 Kalman Filter 估算法更進一步的延伸，

其可將系統未知的參數視為狀態並加以估計，或是處理非線性製程狀態 估計。對於目前濺鍍製程，因沈積率的歷史資料並不完整，故利用 EKF

初始值 P(0)

X

(0) k=0

計算

P ( k + 1 ) k

計算

K k ( + 1)

計算

X ^ˆ ( k + 1 ), ( k Y ^ˆ k + 1 ) k

計算

X ˆ ( k + 1)

計算

P ( k + 1)

k=1,2,…n Q(k)

R(k+1)

Y(k+1)

來估計沈積率的狀態，並將系統參數加入狀態估計。一個狀態空間隨機 離散非線性系統可表示為[23]：

(k + =1) f[ ( ), ( ), ] k w k k ∀ ≥k 0

X X (4.14)

(k + =1) h[ (k +1), (v k+1),k+1] ∀ ≥k 0

Y X (4.15)

X 為系統狀態向量 (n× ，1)

f

與

h

為非線性向量函數，Y 為系統觀測向量

(m× ；w、v 與其他之假設皆與上一節相同。所以(4.14)式與(4.15)式在1) 每一個時期作泰勒展開後，取其一次項之線性化狀態空間方程可寫為：

(k + =1) ( ) ( )k k + ( ) ( )k w k

X F X Γ (4.16)

(k + =1) (k+1) (k + +1) v k( +1)

Y H X (4.17)

其中

( ) ˆ( ), ( ) 0

[ ( ), ( ), ]

( ) ( ) _{X k X k w k}

f k w k k

k k _{= =}

= ∂

∂ F X

X (4.18)

( 1) ˆ( 1 ), ( 1) 0

[ ( 1), ( 1), 1]

( 1)

( 1)

_{X k X k k v k}

h k v k k

k k

_{+ = + + =}

∂ + + +

+ = ∂ +

H X

X

^(4.19)

( ) ˆ( ), ( ) 0

[ ( ), ( ), ]

( ) ( )

_{X k X k w k}

f k w k k

k w k

_{= =}

= ∂

∂

Γ X

(4.20)

則相對應之 Extend Kalman Filter 遞迴方程為：

ˆ(k + =1) f[ ( ), ]ˆ k k + (k +1)[ (k + −1) h[ (ˆ k +1 ),k k+1]]

X X K Y X (4.21)

( k + = 1) ( k + 1 ) k

^T

( k + 1)[ ( k + 1) ( k + 1 ) k

^T

( k + + 1) R k ( + 1)]

⁻1

K P H H P H

(4.22)

( k + 1 ) k = ( ) ( ) k k

^T

( ) k + ( ) ( ) k k

^T

( ) k

P F P F Γ Q Γ

(4.23)

+ = − + + +

ˆ (0) 0, (0) var (0) = =

X P X

(4.25)

其計算流程與 4-1 節相同

4-3 Non-fixed Interval Extended Kalman Filter

事實上，目前量產產品之濺鍍製程中要做到等間距量測是一件困難 的工程，廠商無不想盡辦法提昇工作機台的利用率，為了要使控制器更 能符合實際之情況，就必須把非等間距量測之因素考慮進去。本論文中 是將未量測到之沈積率批次資料，利用預測 d 個時期（d 為未量測之批 次數）之後的沈積率預測值與新量測而得的沈積率資料來更新模型的狀 態及參數，以此類推，求得最佳的沈積率預測值，進而換算出更精確的 沈積時間。以下將說明預測方法及流程。

以 ARIMA(1,1,0)模型為例子，考慮製程中，製程干擾未必屬於 White-Noise，故提高 MA 階次，消除製程干擾之相關性，所以模型轉 變為 ARIMA(1,1,1)，即

( ) (1 ) ( 1) ( 2) ( ) ( 1)

X k = +

φ

X k− −

φ

X k− +a k −

θ

a k− (4.26)

( ) ( ) ( )

Y k = X k + v k

(4.27) 其中X(k)為沈積率狀態，Y(k)為輸出沈積率，a(k)與v(k)分別為濺鍍製程 干擾與量測干擾，ψ與θ為製程模型參數，由(4.26)式與(4.27)式可轉換 為狀態空間方程形式，並將模型參數加入估計。

( ) k = ( k − + 1) a k ( )

X FX Γ

(4.28)

( ) ( ) ( )

Y k = HX k + v k

(4.29) 其中

1

2

( ) ( ) ( )

( ) ( ) X k k X k

k k φ θ

⎡ ⎤

⎢ ⎥

⎢ ⎥

= ⎢ ⎥

⎢ ⎥

⎣ ⎦

X

_，

1 1

1 1 ˆ (0) 0 0 ˆ (0) 0 0 0 1 0 0 0 0 1

X X

φ

φ

⎡ + ⎤

⎢ ⎥

− −

⎢ ⎥

= ⎢ ⎥

⎢ ⎥

⎢ ⎥

⎣ ⎦

F ，

1

0 0 θ

⎡ ⎤

⎢ − ⎥

⎢ ⎥

= ⎢ ⎥

⎢ ⎥

⎣ ⎦

Γ

，

[

^{1 0 0 0}

]

=

H

將上述狀態方程利用卡曼濾波方程預測靶材之沈積率，對於非等間距問 題，處理方法如下：假設兩個非等間距量測點之間(k, k+d)經過了 d 個 批次，利用在第 k 批次時預測第 k+d 批次之沈積率狀態 ˆ (X k +d k)與狀 態變異

P ( k + d k )

並計算在(k+d)批次預測之沈積率 ˆ (Y k +d k)，待新的 沈積率量測值得到後，更新沈積率狀態，計算流程如圖 5 所示，其餘計 算方式皆與 4-2 節相同。

圖 5 Non-Fixed Interval EKF 計算流程 初始值 P(0)

X

(0) k=0

計算

P ( k + 1 ) k

計算

K ( k + 1)

計算

X ˆ ( k + d )

計算

P ( k + d )

k=1,2,…n

Q(k)

R(k+1)

Y(k+d)

計算

X ^ˆ ( k + d k ), ( P k + d k ), ( Y ^ˆ k + d k )

d 為未量測之 lot 數

預測 d 個 lot 以後之 狀態

X 與狀態誤差 P

第五章 建立製程控制與補償機制

在本章中，第一節將介紹沈積率歷史資料的分析方法以及製程中沈積率 變化趨勢，再利用時間序列分析法建立濺鍍製程模型，其次第一節將介紹 如何建立製程模型來預測靶材沈積率以及濺鍍製程之模擬控制。第二、第 三節將介紹模擬與控制之流程與其結果。

5-1 分析濺鍍製程資料

本論文目的在於金屬濺鍍製程沈積率之批次控制，目前是使用半導 體廠的金屬濺鍍製程資料作分析。製程資料方面是由半導體廠鋁製程靶 材的沈積率，皆由同一厚度製程(Al, 4000Å)、同一機臺(Applied Material ENDURA /ENDURA 2)、同一腔體(Chamber 2 or Chamber 3)以及固定輸 入功率之下所量測而得。此濺鍍機台共有4 個製程腔體，Chamber 2 與 Chamber 3 為鋁(Al)濺鍍製程，Chamber C 與 Chamber D 為氮化鈦(TiN) 濺鍍製程，鋁製程之沈積厚度目標值為4000Å，氮化鈦製程之目標厚度 為650Å。本文選用鋁製程資料來討論之原因有二，

1. 鋁製程沈積厚度較氮化鈦厚，靶材消耗量快，故在收集沈積率 資料時，完整的靶材沈積率資料較多。

2. 經由工程師之建議，選用製程較穩定之鋁製程來建立批次控制 機制。

沈積率是由量測機台量測沈積金屬厚度之後，除以製程之沈積時間 而得，在沈積率資料中兩個腔體各有4 個靶材之數據，沈積率歷史資料 如圖6 所示。圖中縱軸為靶材沈積率，橫軸為靶材壽命(Target Life)，一 個靶材約使用 1600 KWH 後就會更換，而且隨著靶材使用時間的增 加，沈積的速率漸漸減少；在同一腔體之其他靶材皆有此現象。因廠商 在製造時為非等間距量測，所以在圖中沈積率資料點的間距呈現大小不 一的現象。為了要分析製程特性，故假設沈積率數列為等間距後再利用 時間序列來擬合。

400 600 800 1000 1200 1400 1600 1800

160 165 170 175

Depo. Rate (Angstrom/sec)

Target Life (KWH) Chamber 2 Target #1

Actual Rate

圖 6 Chamber 2 Target #1 靶材沈積率歷史資料

由圖 6 靶材沈積率隨著使用 KWH 數增加而遞減現象，所以若要建 立濺鍍製程之數學模型，依據第三章之建立時間序列模型方法，先將沈 積率歷史資料先作一次差分後，由 ACF 與 PACF 圖形以及計算各種不 同階次組合之 AIC 值來確定沈積率歷史資料之最佳模式。數據處理過 程如圖7 所示。

圖7 擬合時間數列流程圖

圖中令Y k( ), Y k( −1), Y k( − " 為製程沈積率數列，將其作一次差2) 分，使數列成為平穩數列，差分後的沈積率數列可表示為

( ), ( 1), ( 2) W k W k− W k− "

其中　

( ) ( ) ( 1) (1 ) ( ) W k =Y k −Y k− = −B Y k

沈積率歷史資料 Y(k)

W(k) 一次差分

W k( ) 減去平均值

µ

W

利用時間序列 擬合數列

W k ( )

將差分數列W(k)減去 W(k)之平均值

µ

_W 後，再利用時間數列擬合，即可 得到濺鍍製程時間數列之模型。

Chamber 2 與 Chamber 3 資料擬合後的階次、係數與 AIC 值整理在 表1 與表 2 中。表 1、表 2 中 Target #1 沈積率資料批次數皆為 28 個，

而所擬合的最佳模型階次以及 AIC 值皆不同，可知兩個腔體的製程特 性是不相同的。表 1 中 Chamber 2 的模型階次為 ARIMA(1,1,0)與 ARIMA(2,1,0) ， 而 Chamber 3 的 模 型 階 次 為 ARIMA(2,1,0) 與 ARIMA(4,1,0)，所以在預測靶材沈積率方面選用的時間數列模型為 ARIMA(3,1,0)。

決定此階次的理由：因為各靶材資料批次數不相同，且為非等間距 數列，Target #3~#4 的沈積率資料批次數較少，無法較完整解釋製程的 現象，在AR 項 4 階與 3 階之 AIC 值相近；所以考慮靶材沈積率批次資 料數目以及ARIMA 模型高階次可涵蓋低階次特性，選用 ARIMA(3,1,0) 模型當作是濺鍍製程模型。

表1 Chamber 2 靶材階次、係數與 AIC 值

Target #1 Φ1 Φ2 Φ3 Φ4 θ AIC

(1,1,0) -0.2094 55.79569

(2,1,0) -0.28947 -0.22595 56.46499

(3,1,0) -0.32472 -0.27825 -0.11767 58.12489 (4,1,0) -0.36087 -0.32841 -0.20715 -0.1891 59.16742 沈積率資料批次數：28

Target #2 Φ1 Φ2 Φ3 Φ4 θ AIC

(1,1,0) -0.20213 85.61453

(2,1,0) -0.30914 -0.53891 76.35225 (3,1,0) -0.36575 -0.57194 -0.11209 77.90869 (4,1,0) -0.36219 -0.55748 -0.10216 -0.0263 79.88442 (2,1,1) -0.18249 -0.514 0.17678 77.96719 沈積率資料批次數：37

表2 Chamber 3 靶材階次、係數與 AIC 值

Target #1 Φ1 Φ2 Φ3 Φ4 θ AIC

(1,1,0) -0.39403 50.97645

(2,1,0) -0.59692 -0.40476 47.81872 (3,1,0) -0.68892 -0.51308 -0.1734 48.95098 (4,1,0) -0.74822 -0.68572 -0.35854 -0.25432 49.0519 沈積率資料批次數：28

Target #2 Φ1 Φ2 Φ3 Φ4 θ AIC

(1,1,0) -0.07273 84.47452

(2,1,0) -0.08628 -0.11925 86.05284 (3,1,0) -0.10799 -0.0966 -0.39012 81.08449 (4,1,0) -0.25249 -0.12785 -0.36665 -0.46654 80.17065 (1,1,1) 0.46949 0.92925 82.00692 (2,1,1) 0.63959 -0.28899 0.90562 82.61197 沈積率資料批次數：36

事實上，經由時間序列分析所得到的模型，即為濺鍍製程中的干擾 模型。濺鍍製程控制在目前產業上可分為兩種類型，一種是調整功率，

固定沈積時間來控制沈積率，另外一種則是調整沈積時間而固定輸入功 率來達到產品厚度的穩定。在下一節中將會討論兩種類型的控制方法。

5-2 預測

在半導體方面，已有一些學者[9][10]運用 EWMA 以及 PCC 控制器 來預測濺鍍製程中靶材的沈積率，進而控制沈積時間，使輸出厚度在目 標範圍內。在處理製程漂移現象上DEWMA 已被證明比 EWMA 或 PCC 控制器之效能來的好[24]，因此將 DEWMA 與卡曼濾波器兩者做一比 較，如下表3。其中若製程模型為非線性時變方程，則 Kalman Filter 仍 可適用；在權重的部分，Kalman Gain 功能類似於 DEWMA 之權重，但 Kalman Gain 能夠隨製程狀態調整大小，DEWMA 則否。另一個最大特 點為Kalman Filter 是以狀態空間方程式衍生出之狀態估計方法，所以將 製程干擾與量測干擾都加入考慮，而 DEWMA 只考慮量測干擾，與實 際製程情況有所出入。

表 3 DEWMA 與 Kalman Filter 之比較

估計器類型 DEWMA Kalman Filter

適用範圍 Linear Equation Linear or Non-Linear

參數類型 系統參數為固定值 系統參數可變動

Gain 權值(Weighting)為固定 隨製程而改變

誤差項 只有一個量測誤差項 考慮系統干擾以及量測誤差

參數估計 無 有

基於表 3 之結果，再加上靶材在製程中會做更換，此時 DEWMA 之固

定權重對於下一個靶材不一定適用，故在本篇論文中選用Kalman Filter 來做沈積率的估計。所以濺鍍製程模型假設為：

( ) ( ) ( )

Y k = +

α β

u k +

ε

k (5.1)

其中Y k 為第 k 批次沈積率輸出值， ( )( ) u k 為控制輸入功率，

α

為截距，

β 為製程轉換係數，

ε

( )k 為系統干擾，依據ARIMA(3,1,0)模型，

ε

( )k 又 可寫為

2 3

1 2 3

(1−B)(1−

φ

B−

φ

B −

φ

B ) ( )

ε

k =a k( ) (5.2)

所以將(5.1)式改寫成

2 3

1 2 3

( ) ( ) 1 ( )

(1 )(1 )

Y k u k a k

B B B B

α β φ φ φ

= + +

− − − − (5.3)

現在先對製程干擾項作處理，在半導體業是調整沈積時間來控制產品厚 度，也就是控制器只要能準確的預測下一個批次的沈積率，就能求得下 一個批次的輸入沈積時間。因此在固定輸入功率情況下將(5.3)式控制輸 入項X k 視為常數，則系統模型變為： ( )

2 3

1 2 3

( ) ( ) 1 ( )

(1 )(1 )

Y k k a k

B B B B

ε φ φ φ

= =

− − − − (5.4) 其中

( ) ( ) ( )

Y k =Y k − −

α β

u k (5.5)

(5.4)式表示每批次沈積率的改變是製程干擾所造成，其中Y k 表示沈積( ) 率減去截距以及控制輸入項之後的數列，所以將(5.4)式轉換成狀態空間 形式

1

2 1

3 2

3

1 1 0 0 1

0 1 0 0

( ) ( 1) ( 1)

0 0 1 0

0 0 0 0

k k a k

φ φ φ φ φ

φ

⎡ + ⎤ ⎡ ⎤

⎢ − ⎥ ⎢ ⎥

⎢ ⎥ ⎢ ⎥

= − + −

⎢ − ⎥ ⎢ ⎥

⎢ − ⎥ ⎢ ⎥ ⎣ ⎦

⎣ ⎦

X X

(5.6)

[ ]

( ) 1 0 0 0 ( ) ( )

Y k =

X

k + v k

(5.7)

式中X

( ) k

為第 k 批次沈積率之狀態，

Y k ( )

為第 k 批次處理過後之輸出 沈 積 率 ，

a k ( − 1)

、

v k ( )

分 別 為 製 程 與 量 測 干 擾 ， 兩 者 皆 屬 於 White-Noise。再將 ARIMA 模型係數當作被估計狀態放入狀態方程式 中，以4-2 節線性化方法可得到線性化之狀態空間方程

( ) ( 1) ( 1) ( ) ( ) ( )

k k a k

Y k k v k

= − + −

= +

X FX Γ

HX (5.8)

其中

1

2

3

4

1

2

3

( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) X k X k X k k X k k k k

φ φ φ

⎡ ⎤

⎢ ⎥

⎢ ⎥

⎢ ⎥

⎢ ⎥

= ⎢ ⎥

⎢ ⎥

⎢ ⎥

⎢ ⎥

⎢ ⎥

⎣ ⎦

X ，

1 1

2 1 1 1

3 2 1 1

3 1

1 1 0 0 ˆ (0) 0 0

ˆ ˆ

0 1 0 (0) (0) 0

ˆ ˆ

0 0 1 0 (0) (0) 0 0 0 0 0 ˆ (0)

0 0 0 0 1 0 0

0 0 0 0 0 1 0

0 0 0 0 0 0 1

X

X X

X X

X

φ

φ φ φ φ

φ

⎡ + ⎤

⎢ ⎥

− −

⎢ ⎥

⎢ ⎥

− −

⎢ ⎥

= ⎢⎢ − − ⎥⎥

⎢ ⎥

⎢ ⎥

⎢ ⎥

⎢ ⎥

⎣ ⎦

F ，

1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

⎡ ⎤

⎢ ⎥

⎢ ⎥

⎢ ⎥

⎢ ⎥

= ⎢ ⎥

⎢ ⎥

⎢ ⎥

⎢ ⎥

⎢ ⎥

⎣ ⎦

Γ ，H=

[

1 0 0 0 0 0 0

]

，

利用(5.8)式與 Extend Kalman Filter 方法，即可計算下一個批次的沈 積率，進而計算出下一個批次的沈積時間。沈積時間計算公式為：

Target

Predicted

Thickness Deposition Time

Deposition Rate

= (5.9)

為了測試Extend Kalman Filter 預測之效能，以 Chamber 2 與 Chamber 3 的沈積率歷史資料為例，將Extend Kalman Filter 與 DEWMA 兩種方法 作比較，其中DEWMA 控制器之預測形式寫為[24]

ˆ( ) ˆ( 1) ˆ( 1) ( )

Y k =a k− + p k − +bu k (5.10) 更新截距與漂移項為

1 1

2 2

ˆ( ) ( ( ) ( )) (1 )( (ˆ 1) ˆ( 1)) ˆ( ) ( ( 1) ( 1) ˆ( 1)) (1 ) (ˆ 1) a k w Y k bu k w a k p k

p k w Y k bu k a k w p k

= − + − − + −

= − − − − − + − − (5.11) 控制法則為

ˆ ˆ Target ( ) ( ) ( 1) a k p k

u k b

− −

+ = (5.12)

其利用估計 k 時期製程之截距與漂移量，將之當作為下一個時期(k+1) 的截距與漂移量，再加上k+1 時期之控制輸入 u(k+1)來預測 k+1 時期之 沈積率。若要求得 k+1 時期之控制輸入值，可利用(5.11)式所估計之截 距與漂移項，以及(5.10)式之預測值

Y k ˆ( )

改為製程目標值求得。DEWMA 之權重設定如表 4 所示，EKF 所用之參數是取自表 1 Target #1 ARIMA(3,1,0)之參數，在效能比較上選用 MSE(Mean Square Error)來作 指標，其計算方式為

2

1

( ( ) ˆ ( ))

n

k

Y k Y k

MSE n

=

−

= ∑

(5.13)

兩種預測方法之 MSE 比較值列於表 2、3 中。沈積率在 Chamber 2 中 預測結果如圖8，Chamber 3 預測結果如圖 9。

表4 DEWMA 權重

DEWMA w1 w2 Chamber 2 0.6 0.6

Chamber 3 0.6 0.6

表 5 Chamber 2 兩種預測方法之 MSE 值 Target #1 Target #2 Target #3 Target #4 DEWMA 1.13725 1.39385 3.36120 2.17500

EKF 1.34397 1.37776 3.10772 2.55840

表 6 Chamber 3 兩種預測方法之 MSE 值 Target #1 Target #2 Target #3 Target #4 DEWMA 1.12579 1.76262 1.82371 3.37897

EKF 0.97490 1.38969 2.83829 3.25083

400 600 800 1000 1200 1400 1600 1800

160 162 164 166 168 170 172 174 176

Target Life (KWH)

Depo. Rate (Angstrom/sec)

EKF & DEWMA Chamber 2 Target #1

Actual Rate EKF Predicted DEWMA Predicted

圖 8 兩種預測方法對 Chamber 2 Target #1 之沈積率預測結果比較圖

0 200 400 600 800 1000 1200 164

166 168 170 172 174 176

Target Life (KWH)

Depo. Rate (Angstrom/sec)

EKF &DEWMA Chamber 3 Target #1

Actual Rate EKF Predicted DEWMA Predicted

圖 9 兩種預測方法對 Chamber 3 Target #1 之沈積率預測結果比較圖 由以上圖形及 MSE 值可以得知，實際沈積率數據並不完整，所以造成 在Target #3~#4 預測沈積率時，發生預測與實際沈積率之 MSE 值較大 現象。同時也表示ARIMA(3,1,0)模型不適用於 Target #3~#4，當沈積率 批次資料越多時，EKF 會比 DEWMA 更準確預測沈積率。

為解決非等間距問題，利用 4-4 節之方法加以修改 Extend Kalman Filter，將未量測的批次數考慮進去進行預測，另與 A. Chen, and R. S.Guo, 所提出的Age-Based DEWMA Controller[14]比較，Age-Based DEWMA 之預測與控制方法與 DEWMA 類似，但其加入兩個兩側點中間未量測 的批次數來修正漂移量估計與控制輸入值，其預測沈積率形式與(5.10) 式相同

ˆ( ) ˆ( 1) ˆ( 1) ( )

Y k =a k− + p k − +bu k (5.10) 估計截距與漂移量表示為

1 1

2 2

ˆ( ) ( ( ) ( )) (1 )( (ˆ 1) ( ) (ˆ 1)) ( ( 1) ( 1) ˆ( 1))

ˆ( ) (1 ) (ˆ 1)

( )

a k w Y k bu k w a k d k p k Y k bu k a k

p k w w p k

d k

= − + − − + −

− − − − −

= + − − (5.14)

控制律表示為

ˆ ˆ

Target [ ( ) ( 1) ( )]

( 1) a k d k p k

u k b

− + +

+ = (5.15) d(k)為 k 時期與 k-1 時期量測點中間，未量測之批次數，Age-Based DEWMA Controller 之權重是利用 Chamber 2 Target #1 的沈積率資料計 算出最小的 MSE 而得。權重值列於表 7，預測結果如圖 10、11，兩種 預測方法之 MSE 值列於表 8、9 中。將等間距與非等間距預測之 MSE 值列於表10、11。

表7 Age-Based DEWMA 之權重值 Age-Based DEWMA w1 w2

Chamber 2 0.69 0.98 Chamber 3 0.69 0.98 表8 Chamber 2 兩種預測之 MSE 值

MSE Target #1 Target #2 Target #3 Target #4 Age-Based DEWMA 1.36720 2.01495 5.71981 3.17630 Non-Fixed Interval EKF 0.88636 1.08991 3.06715 2.11061

表9 Chamber 3 兩種預測之 MSE 值

MSE Target #1 Target #2 Target #3 Target #4 Age-Based DEWMA 1.81482 1.59617 2.63959 5.03984 Non-Fixed Interval EKF 0.89655 1.38969 1.82715 2.47128

表10 Chamber 2 等間距與非等間距預測結果比較

MSE Target #1 Target #2 Target #3 Target #4 DEWMA 1.13725 1.39385 3.36120 2.17500 Age-Based DEWMA 1.36720 2.01495 5.71981 3.17630 EKF 1.34397 1.37776 3.10772 2.55840 Non-Fixed Interval EKF 0.88636 1.08991 3.06715 2.11061

表11 Chamber 3 等間距與非等間距預測結果比較

MSE Target #1 Target #2 Target #3 Target #4 DEWMA 1.12579 1.76262 1.82371 3.37897 Age-Based DEWMA 1.81482 1.59617 2.63959 5.03984 EKF 0.97490 1.38969 2.83829 3.25083 Non-Fixed Interval EKF 0.89655 1.38969 1.82715 2.47128

400 600 800 1000 1200 1400 1600 1800

160 162 164 166 168 170 172 174 176

Target Life (KWH)

Depo. Rate (Angstrom/sec)

NFIEKF & Age-Based DEWMA Chamber 2 Target #1

Actual Rate

Non-Fixed Interval EKF Predicted Age-Based DEWMA Predicted

圖 10 兩種預測方法對 Chamber 2 Target #1 之沈積率預測結果比較圖

0 200 400 600 800 1000 1200 162

164 166 168 170 172 174 176 178

Target Life (KWH)

Depo. Rate (Angstrom/sec)

NFIEKF & Age-Based DEWMA Chamber 3 Target #1

Actual Rate

Non-Fixed Interval EKF Predicted Age-Based DEWMA Predicted

圖11 兩種預測方法對 Chamber 3 Target #1 之沈積率預測結果比較圖 分析兩個Chamber 的預測結果圖與 MSE 值可以歸納出以下幾點：

1. 經過非等間距方法修正過後的 EKF，其預測效能較其他預測方法優 良。

2. 在 Target #3~4 的沈積率預測值都與實際沈積率相差甚大，此因沈積 率資料流失過多，以致於製程模型與實際製程產生誤差。

3. Age-Based DEWMA 控制器之權重 w2 趨近於 1，表示此控制器是以 近乎Random Walk 的型態來進行預測。

4. 由 MSE 表格可知，若使用目前業界常用的 DEWMA 控制器，其權 重會隨著更換靶材或其他因素而改變，而EKF 方法可克服此一缺點。

5-3 製程模擬與控制

因為半導體廠目前以產能為第一優先，無法實際上機驗證控制方 法，所以利用模擬製程來驗證控制方法是否合適，並與 DEWMA 控制 器比較效能。

5-3-1 模擬方法

濺鍍製程的模擬是利用(5.4)式作為製程模型，以 5-1 節所決定之 ARIMA(3,1,0)型態加入製程模型，並固定輸入功率。則模型根據 5-1 節 將數據差分一次後，可寫為

2 3

1 2 3

(1 − φ B − φ B − φ B )( ( )) W k = a k ( )

(5.16) 其中

( ) ( )

_w

W k

=

W k

−

µ

(5.17)

[ ( )]

_W

E W k = µ

(5.18)

在模擬的步驟中，所用的數值皆由濺鍍製程之歷史資料 Chamber 3 Target #2 計算而得，模擬步驟整理如下：

1. 設沈積率目標值為 176(Å/sec)，產生 300 個批次沈積率數據

2. 取用 Chamber 3 Target #2 ARIMA(3,1,0)係數

φ

₁、

φ

₂、

φ

₃，作為模擬模 型的係數。

3. 給定

a k ( ) ~ (0, N σ

_a²

)

，經由(5.19)式產生 W(k)數列，給定起始沈積率

後再還原成模擬沈積率

Y k ( )

2 3

1 2 3

( ) ( )

(1 )

^w

W k a k

B B B µ

φ φ φ

= +

− − − (5.19) 經由模擬產生的沈積率數據如圖12。

0 500 1000 1500

158 160 162 164 166 168 170 172 174 176 178

Simulation & Chmaber 3 Target 2 Data

Target Life (KWH)

Depo. Rate (Angstrom/sec)

Simulation Actul Rate

圖12 模擬產生沈積率與實際沈積率比較圖

5-3-2 控制方法

光碟片濺鍍製程是調整製程的輸入功率而固定沈積時間來控制沈 積率以及薄膜厚度的穩定度，控制架構如圖13 所示

圖13 濺鍍控制架構圖

圖中 u(k)為輸入功率，Y(k)為沈積率，e(k)為輸出沈積率與目標值之誤 差，T 為沈積率目標值，

ε ( ) k

為製程干擾；圖中上方濺鍍製程是由 5-1 節所提出的模擬製程模型，中間之 EKF 主要功能為預測下一個批次的 干擾以及估計製程參數，提供給控制器計算下一個批次所需設定的功率 值，下方控制器主要利用 EKF 以及目標值之資訊計算出下一個輸入功 率值。濺鍍製程模擬與 EKF 估計方法在前面章節已有說明，接下來要 說明的是控制器部分。

假設製程模型和製程干擾

ε ( ) k

設定與(5.1)式相同

( ) ( ) ( )

Y k

= +

α β u k

+

ε k

(5.1) 為使沈積率在製程中維持在一定準位上，就必須調整控制輸入功率來消

Sputter Process

( ) k

ε

Extended Kalman Filter

( ) ( ) ( ) u k k Y k α β + + ε = ( )

u k

Controller T

u k+ ( 1)

1 2 3

ˆ ˆ ˆ

ˆ( 1| ), ( 1), ( 1), ( 1) k k k k k

ε

+

φ

+

φ

+

φ

+