Daftar Pustaka__________________________________________________________________
DAFTAR PUSTAKA
[1] Farina, A., dan Studer, F. A. (1985), Radar Data Processing, Vol. I –
Introduction and Tracking, Research Studies Press. [2] Farina, A., dan Studer, F. A. (1985), Radar Data Processing, Vol. II –
Advanced Topics and Applications, Research Studies Press. [3] Kayton, M., dan Fried, W. R. (1997), Avionics Navigation Systems, Second
Edition, John Wiley & Sons. [4] R.E. Kalman, R.S Bucy, New results in Linear Filtering and Prediction
Theory, Baltimore, Maryland (1961). [5] Brookner Eli. (1998), Tracking and Kalman Filtering Made Easy, Wiley & Sons, inc, ISBN 0-471-18407. [6] Gustafsson, F., Ljung, L., Millnert, M. (2001), Signalbehandling, Studentlitteratur, pp. 289-336, ISBN 91-44-01709-X. [7] Bar-Shalom, Y., Rong Li, X., Kirubarajan, T., (2001), Estimation with
applications to tracking and navigation. Theory, Algorithms and Software, John Wiley & Sons. [8] Park, S-T., dan Lee J. G. (1998), Design of a Practical Tracking Algoritm with Radar Measurements,
IEEE Transactions on Aerospace and
Electronic Systems Vol. 34, No. 4 (October 1998), halaman 1337 - 1344. [9] Chong, C. Y., Mori, S., Barker, W. H., dan Chang, K. C. (2000), Architectures and Algorithms for Track Association and Fusion, IEEE AES
Systems Magazine (January 2000), halaman 5 – 13. [10]
Bar-Shalom, Y.
dan Fortman, T.E., (1988), Tracking and Data
Association, San Diego, CA Academic Press. [11]
Noren, B. O. (2004), STCA – an aicraft confilict alert system, Master
Thesis, Linkoping Universiteit .
73
Daftar Pustaka__________________________________________________________________
[12]
Welch, G. and Bishop, G. An Introduction to the Kalman Filter, TR
95-041 Department of Computer Science. [13]
Muhammad, H. (2002) HandOut Kuliah PN-553 Identifikasi
Parameter Pesawat Udara,
Departemen Teknik Penerbangan, Institut
Teknologi Bandung. [14]
Frank, O. (February 2003) Multiple Target Tracking, Master Thesis,
Swiss Federal Institute of Technology Zurich (ETH).
74
Lampiran A____________________________________________________________________
LAMPIRAN A Konsep Dasar Estimasi Berikut ini akan diberikan sebuah pengenalan kepada masalah estimasi dasar. Masalah estimasi parameter random (acak) dan non-random akan dijelaskan dengan tambahan beberapa estimator umum.
Masalah estimasi sebuah parameter time invariant dijelaskan sebagai berikut
[14]. Dari observasi z j = h(x j , w j , j ) yang dibuat dengan noise w j , harus
ˆ (z , k ) yang meng-estimasi ditemukan sebuah fungsi dari observasi ˆx k j = X i: j ˆ (.,.) disebut aturan estimasi atau estimator dan nilai dari x pada waktu k. X z i: j ={z i ,...., z j } adalah semua observasi dari waktu i ke waktu j . Estimator Δ
harus didefinisikan dengan memperhatikan suatu kriteria kinerja yang menmberikan perbandingan dari aturan-aturan berbeda dan estimasi yang mungkin. Sekarang akan diperkenalkan suatu fungsi yang disebut fungsi harga
(cost function) κ (ˆx − x ) yang memberikan konsekuensi bila terjadi perbedaan antara nilai sebenarnya dan estimasi. Sering fungsi harga
κ (ˆx − x ) = ˆx − x
2
(A.1)
digunakan.
Perbedaan dari berbagai macam estimator dapat dikelompokkan dalam tiga tipe umum berikut : [7]
75
Lampiran A____________________________________________________________________
• Predictor, menggunakan observasi hanya pada nilai sebelumnya sampai tahap waktu k, dengan tujuan bahwa variabel keadaan dari sistem dinamik akan di-estimasi : i ≤ j < k
• Filters, menggunakan observasi sampai dan termasuk di dalam tahap waktu k, dengan tujuan bahwa variabel keadaan dari sistem dinamik akan diestimasi : i ≤ j = k
• Smoothers, menggunakan observasi lebih dari tahap waktu k, dengan tujuan bahwa variabel keadaan dari sistem dinamik akan di-estimasi : k < i ≤ j
Selanjutnya, estimator dari tipe kedua di atas (Filters) akan dijelaskan, dimana diambil asumsi bahwa i = 1 . Sehingga, masalahnya adalah menemukan aturan estimasi untuk : ˆ (z , k ) xˆ k k = X 1:k
(A.2)
Pada dasarnya, ada dua cara berbeda untuk meng-estimasi sebuah parameter, tergantung dari apakah parameter itu acak (random) atau tidak (non random) :
• Pendekatan Bayesian (parameter acak) : Dimulai dengan a priori probability density function (pdf) (fungsi kepadatan probabilitas) p (x ) dari parameter acak yang a posteriori pdf nya dapat diperoleh dengan menggunakan formula Bayes :
p (x z 1:k ) =
p(z 1:k x ) p( x )
(A.3)
p(z 1:k )
Dua metode yang umum digunakan adalah estimator maximum a posterior (MAP) dan estimator minimum mean-squared error (MMSE).
• Pendekatan Non-Bayesian (parameter tidak acak) :
76
Lampiran A____________________________________________________________________
Berlawanan dari penjelasan di atas, tidak ada a priori pdf yang diasosiasikan dengan parameter. Pada kasus ini, pdf pengukuran dikondisikan sesuai dengan parameter, yang disebut fungsi likelihood λ( x ) dari parameter sebagai suatu ukuran seberapa mirip sebuah nilai parameter, bila diberi observasi yang telah dimiliki : λ(x ) = p (z 1:k x )
(A.4)
Estimator maximum likelihood (ML) dan metode least squares (LS) termasuk pada kelas estimator ini.
C.1. Maximum Likelihood Estimator Aturan yang banyak dipakai untuk estimasi parameter tidak acak adalah metode ML. Tujuannya adalah menentukan nilai dari parameter x dengan memaksimalkan fungsi likelihood λ(x ) . Estimasi ML dijelaskan dengan : ˆx ML = arg max λ(x ) = arg max p (x z 1:k ) x
x
(A.5)
Harus dingat, ketika x adalah konstanta yang tidak diketahui, ˆx ML , sebagai fungsi dari beberapa set observasi acak z 1:k , adalah variabel acak.
C.2. Maximum A Posteriori Estimator Di dalam kasus dimana parameter x dianggap acak, realisasi dari x menurut suatu pdf p ( x ) diasumsikan telah ada. Aturan Bayes digunakan untuk mencari a posterior pdf :
p (x z 1:k ) =
p(z 1:k x ) p( x )
(A.6)
p(z 1:k )
A posteriori pdf p (x z 1:k ) adalah suatu ukuran seberapa mirip sebuah nilai parameter, bila diberi observasi hingga tahap waktu k dan informasi a priori p (x )
Sehingga untuk suatu parameter acak. Estimator-MAP adalah ekivalen dengan estimator-ML untuk parameter tidak acak :
77
Lampiran A____________________________________________________________________
ˆx MAP = arg max p (x z 1:k )
(A.7)
x
= arg max[ p (z 1:k x ) p (x )]
(A.8)
x
Telah jelas, jika a priori pdf p (x ) adalah suatu distribusi yang seragam (disebut juga diffuse prior), maka estimasi ML dan estimasi MAP akan berharga sama. Seperti yang dijelaskan pada [7], “Estimator ML non Bayesian bukan lain adalah estimator MAP Bayesian dengan ketidak-tahuan
a priori yang lengkap, direfleksikan dengan distribusi a priori yang seragam. Pernyataan ini menghasilkan pandangan filosofis yang menyatukan pendekatan Bayesian dan non Bayesian terhadap estimasi. Harus dicatat bahwa estimator-MAP adalah kasus khusus dari estimasi Bayesian dimana nilai yang diharapkan dari fungsi harga
⎧0,ˆx = x κ(ˆx − x ) = ⎨ ⎩1, ˆx ≠ x
(A.9)
Harus diminimalkan. Nilai harapan dari fungsi harga adalah sama dengan probabilitas kesalahan yang diberi observasi z 1:k : E [κ(ˆx − x ) z 1:k ] = P (X ≠ x Z1:k = z 1:k )
(A.10)
Sehingga estimator-MAP didapat : ˆx MAP = arg min P (X ≠ x Z1:k = z 1:k ) x
= arg max P (X = x Z1:k = z 1:k ) x
= arg max P(X = x ) x
p(z 1:k x ) p(z 1:k )
= arg max P(X = x ) p (z 1:k x ) x
(A.11) (A.12) (A.13) (A.14)
C.3. Estimator Least Squares Satu lagi estimator untuk parameter tidak acak adalah metode LS. Diasumsikan bahwa serangkaian pengukuran z 1:k diketahui mempunyai hubungan lewat fungsi penentu h dalam bentuk :
78
Lampiran A____________________________________________________________________
z j = h(x , j ) + w k
(A.15)
Estimasi kuadrat terkecil dari x setelah mengambil k pengukuran adalah nilai yang meminimalkan jumlah dari kesalahan kuadrat : x
[
][
]
ˆx LS = arg min ∑ z j − h(x , j ) T z j − h(x , j ) x
(A.16)
j =1
Tidak ada asumsi yang dibuat untuk noise w k .
C.4. Minimum Mean-Squared Error Estimator Metode yang ekivalen dengan estimator-LS untuk parameter acak didapat dengan meminimalkan kesalahan kuadrat nilai tengah, estimator-MMSE.
[
ˆx MMSE = arg min E (ˆx − x )T (ˆx − x ) z 1:k x
]
(A.17)
Dimana w k adalah serangkaian variabel acak dengan pdf yang diketahui. Solusi dari masalah ini didapatkan dengan menurunkan nilai yang diharapkan terhadap estimasi ˆx :
[
∂ T E (ˆx − x ) (ˆx − x ) z 1:k ∂ˆx =
]
∂ ⎡ ∞ (ˆx − x )T (ˆx − x ) p(x z 1:k )dx⎤⎥ ∫ ⎢ − ∞ ⎣ ⎦ ˆ ∂x
(A.18)
(ˆx − x ) p(x z 1:k )dx = 0
(A.19)
ˆx MMSE = ∫ xp (x z 1:k )dx = E [x z 1:k ]
(A.20)
= 2∫
∞
−∞
Didapatkan : ∞
−∞
Yang berarti, bahwa estimasi MMSE adalah sama dengan nilai tengah bersyarat x yang diberi semua observasi z 1:k . Perhatikan bahwa proses estimasi MMSE adalah suatu kasus khusus dari proses estimasi Bayesian, dimana nilai yang diharapkan dari fungsi harga κ(ˆx − x ) = ˆx − x
2
harus
diminimalkan :
[
ˆx MMSE = arg min E [κ(ˆx − x ) z 1:k ] = E ˆx − x z 1:k x
79
2
]
(A.21)
Lampiran A____________________________________________________________________
Metode-metode di atas dipaparkan untuk memberikan dasar-dasar dan gambaran mengenai masalah estimasi dan estimator sederhana yang sering digunakan.
80
Lampiran B____________________________________________________________________
LAMPIRAN B Penurunan Persamaan Filter Kalman Seperti pada filter Wiener, Filter Kalman dipilih untuk meminimalkan error kuadrat nilai tengah[14]
[(
ˆx i j = arg min n E x i − ˆx i ˆx (i j )∈R
j
)(x
i
− ˆx i
)
T
j
z 1 ,..., z j
]
(B.1)
Seperti telah diturunkan pada Lampiran A, kesalahan kuadrat nilai tengah adalah nilai kondisi nilai tengah dari keadaan x i pada tahap waktu i, berdasarkan semua observasi z 1 ,..., z j yang tergantung dari tahap waktu j.
[
ˆx i j = E x i z 1: j
]
(B.2)
Dimana z 1:k menunjukkan semua observasi sampai tahap waktu k. Estimasi a
priori ˆx k k −1 didapatkan dengan memakai estimasi sebelumnya ˆx k −1 k −1 seperti dalam persamaan (2.8) pada Bab II:
ˆx k k −1 = E [x k z 1:k −1 ]
(B.3)
= E [Fk x k −1 + B k u k + G k v k z 1:k −1 ]
(B.4)
= Fk E [x k −1 z 1:k −1 ] + B k u k + G k E [v k z 1:k −1 ]
(B.5)
= Fk ˆx k −1 k −1 + B k u k
(B.6)
Estimator linier yang rekursif ingin dicari untuk keadaan x(⋅) . Jadi proses estimasi dapat ditulis sebagai jumlah linier yang diberi bobot dari prediksi dan observasi yang baru :
81
Lampiran B____________________________________________________________________
ˆx k k = K k ˆx k k −1 + K k z k
(B.7)
Dimana matriks K k dan K k belum diketahui. Masalah optimasi ini dapat dipecahkan dengan memakai prinsip orthogonalitas. Kondisinya ditulis sebagai :
[( E [(x
) ] )z ] = 0
E x k − ˆx k k z Ti = 0 k
− ˆx k k
i = 1,..., k − 1
(B.8)
T k
(B.9)
Dengan men-substitusi x k dari (2.8) dan ˆx k k dari (B.7) maka didapatkan
[(
) ]
E Fk x k −1 + B k u k + G k v k − K k ˆx k k −1 + K k z k z Ti = 0
(B.10)
Dan dengan mengganti z k sesuai (2.9)
[(
) ]
E Fk x k −1 + B k u k + G k v k − K k ˆx k k −1 + K k H k x k + K k w k z Ti = 0 (B.11)
[
]
Fk E [x k −1 ]z Ti + −K k E ˆx k k −1 z Ti − K k H k E [x k ]z Ti + B k u k = 0 (B.12) Dimana fakta bahwa E [w k ] = E [v k ] = 0 telah digunakan. Kemudian dengan men-substitusi x k lagi dari (2.8) maka
[
]
Fk E [x k −1 ]z Ti − K k E ˆx k k −1 z Ti − K k Fk H k E [x k ]z Ti = 0 (B.13)
Fk x k −1 dapat diekspresikan sebagai x k − B k u k − G k v k dengan menggunakan (2.8) dan didapatkan
82
Lampiran B____________________________________________________________________
[
]
E [x k ]z Ti − K k E ˆx k k −1 z Ti − K k H k E [x k ]z Ti = 0
[
(B.14)
)]
(
E (x k − K k H k x k − K k x k ) − K k ˆx k k −1 − x k z Ti = 0
(B.15)
Dengan menggunakan (B.6) dan (2.8) maka
[
]
[
E x k − ˆx k k −1 z Ti
]
= E Fk x k −1 + G k v k − Fk x k k −1 z Ti
(B.16)
= Fk E x k − ˆx k −1 k −1 z Ti
(B.17)
[(
) ]
Dan bila dibandingkan dengan (B.9) maka dapat dilihat bahwa pernyataan ini adalah nol. Jadi akhirnya didapat
(I − K
k
H k − K k )E [x k ]z Ti = 0
(B.18)
Persamaan ini dapat dipakai oleh setiap x k jika
Kk = I − KkHk
(B.19)
Disubstitusikan ke dalam (B.7) menghasilkan
ˆx k k = (I − K k H k )ˆx k k −1 + K k z k
(
(B.20)
= ˆx k k −1 + K k z k − H k ˆx k k −1
)
(B.21)
Harus diperhatikan bahwa H k ˆx k k −1 adalah observasi yang diprediksi, sehingga bentuk terakhir dari persamaan estimasi ini dapat diinterpretasikan sebagai penjumlahan dari prediksi ditambah bagian K k dari perbedaan antara
83
Lampiran B____________________________________________________________________
yang sebenarnya dengan observasi yang diprediksi. Gain K k masih harus ditentukan.
Kovarian kesalahan estimasi didefinisikan sebagai kesalahan kuadrat nilai tengah dari estimasi :
Δ
[(
Pi j = E x i − ˆx i
j
)(x
i
− ˆx i
)
T
j
z 1: j
]
(B.22)
Didapat matriks kovarian kesalahan estimasi a priori dengan mengurangkan (B.6) dengan (2.8) :
(
)
ˆx k − ˆx k k −1 = Fk x k −1 − ˆx k −1 k −1 + G k v k
(B.23)
Dan mengambil nilai harapan yang telah dikondisikan pada semua observasi sampai pada waktu k-1, menghasilkan
[(
)(
)
T
Pk k −1 = E x k − ˆx k k −1 x k − ˆx k k −1
[(
z 1:k −1
)(
= Fk E x k − ˆx k k −1 x k − ˆx k k −1
)
T
]
(B.24)
]
z 1:k −1 FkT + G k Q k G Tk (B.25)
= Fk Pk −1 k −1FkT + G k Q k G Tk
(B.26)
Lebih jauh, matriks kovarian kesalahan estimasi a posteriori didefinisikan
[(
)(
Pk k = E x k − ˆx k k −1 x k − ˆx k k −1
)
T
z1:k
]
(B.27)
Selanjutnya, (2.9) disubstitusikan ke dalam (B.21) dan kemudian mensubstitusikan hasilnya untuk ˆx k k ke dalam (B.27) menghasilkan
84
Lampiran B____________________________________________________________________
[(
(
))(
(
Pk k = E x k − ˆx k k −1 + K k H k ˆx k + w k − H k ˆx k k −1 x k − ˆx k k −1 + K k H k ˆx k + w k − H k ˆx k k −1 (B.28)
Setelah itu dengan melakukan proses ekspektasi dan memperhatikan bahwa
[ E [x
E x k − ˆx k k −1 z 1:k k
− ˆx k k −1 w Tk
] adalah kovarian ] = 0 , didapatkan
kesalahan estimasi a priori Pk k −1 dan
Pk k = (I − K k H k )Pk k −1 (I − K k H k ) + K k R k K Tk T
(B.29)
Persamaan di atas adalah ekspresi umum dari kovarian kesalahan estimasi a
posteriori dan berlaku untuk semua harga K k , entah itu optimal atau tidak.
Kembali kepada masalah optimasi, diharapkan dapat ditemukan harga K k yang khusus untuk meminimalkan kesalahan estimasi kuadrat nilai tengah yang berkondisi :
[(
)(
E x k − ˆx k k −1 x k − ˆx k k −1
( [( = trace(P )
)
T
)(
z 1:k
]
= trace E x k − ˆx k k −1 x k − ˆx k k −1
)
T
z1:k
])
(B.30) (B.31)
k k
Persyaratan individual pada diagonal mayor dari Pk k harus diminimalkan karena persyaratan ini merepresentasikan variansi kesalahan estimasi untuk elemen-elemen dari vektor keadaan yang sedang diestimasi. Untuk melakukan optimasinya diperlukan dua formula matriks differensial :
∂ [trace(AB )] = B T ∂A
[
]
(A, B ) harus matriks persegi
∂ trace(ACA T ) = 2 ACT ∂A
C harus matriks simetrik
85
(B.32) (B.33)
)) ] T
Lampiran B____________________________________________________________________
Bentuk umum Pk k dari (B.29) harus dikembangkan dan ditulis :
(
)
Pk k = Pk k −1 − K k H k Pk k −1 − Pk k −1 H Tk K Tk + K k H k Pk k −1 H Tk + R k + K Tk (B.34)
Dengan menurunkan trace dari matriks Pk k terhadap K k dan mengaplikasikan dua formula matriks differensial, maka didapatkan
[
]
(
∂ tracePk k −1 = −2 H k Pk k −1 ∂K k
)
T
(
)
+ 2K k H k Pk k −1 H Tk + R k (B.35)
Dengan mengingat fakta bahwa trace dari Pk k −1 H Tk K Tk adalah sama dengan trace dari transpose nya K k H k PkT k −1 = K k H k Pk k −1 , penurunannya kemudian diset berharga nol dan memberikan pemecahan untuk gain optimal K k . Sehingga didapatkan
(
K k = Pk k −1 H Tk H k Pk k −1 H Tk + R k
)
−1
(B.36)
Dengan harga K k yang khusus ini, yang dinamakan Kalman Gain, persamaan (B.21) menjadi estimator kesalahan kuadrat nilai tengah minimal linier yang optimal
Persamaan-persamaan (B.36), (B.21), dan (B.29) yang telah didapat di atas secara berurutan dipakai pada persamaan-persamaan (2.17), (2.18), dan (2.19) di Bab II.
86
Lampiran C____________________________________________________________________
LAMPIRAN C Aplikasi Praktis Filter Kalman Pada Simulasi Radar Tracking C.1. Penurunan Matriks Transisi Keadaan Φ
Bila diambil tipe gerak dengan kecepatan konstan sesuai dengan persamaan (2.4), maka fungsi keadaan f (x ) berisi persamaan-persamaan dari variabel keadaan x , yaitu
⎡ x1 ⎤ ⎡ x ⎤ ⎢x ⎥ ⎢V ⎥ x x = ⎢ 2⎥ = ⎢ ⎥ ⎢x 3 ⎥ ⎢ y ⎥ ⎢ ⎥ ⎢ ⎥ ⎣x 4 ⎦ ⎢⎣V y ⎥⎦
(C.1)
⎡ x k +1 ⎤ ⎡ x k + (V x ,k ⋅ Δt )⎤ ⎡ x + (V x ⋅ Δt )⎤ ⎢V ⎥ ⎢ ⎢ ⎥ ⎥ V x ,k Vx x ,k +1 ⎥ ⎥ Æ f (x ) = ⎢ ⎥ dimana ⎢ =⎢ ⎢ y k +1 ⎥ ⎢ y k + (V y ,k ⋅ Δt )⎥ ⎢ y + (V y ⋅ Δt )⎥ ⎢ ⎥ ⎢ ⎢ ⎥ ⎥ V y ,k Vy ⎣⎢V y ,k +1 ⎦⎥ ⎣⎢ ⎣⎢ ⎦⎥ ⎦⎥
(C.2)
Dari persamaan (C.2), matriks transisi keadaan Φ diturunkan seperti di bawah ini :
⎡ ∂F1 ⎢ ∂x ⎢ 1 ⎢ ∂F1 ∂f (x ) ⎢ ∂x 2 =⎢ Φ= ∂F ∂x ⎢ 1 ⎢ ∂x 3 ⎢ ∂F1 ⎢ ⎣ ∂x 4
∂F2 ∂x1 ∂F2 ∂x 2 ∂F2 ∂x 3 ∂F2 ∂x 4
∂F3 ∂x1 ∂F3 ∂x 2 ∂F3 ∂x 3 ∂F3 ∂x 4
∂F4 ⎤ ∂x1 ⎥ ⎥ ∂F4 ⎥ ∂x 2 ⎥ ∂F4 ⎥ ⎥ ∂x 3 ⎥ ∂F4 ⎥ ⎥ ∂x 4 ⎦
Sehingga didapatkan
87
(C.3)
Lampiran C____________________________________________________________________
⎧1 Δt ⎪0 1 ⎪ Φ=⎨ ⎪0 0 ⎪⎩0 0
0 0⎫ 0 0 ⎪⎪ ⎬ 1 Δt ⎪ 0 1 ⎪⎭
(C.4)
Persamaan (C.4) di atas dipakai pada persamaan (3.2) di Bab III.
C.2. Penurunan Matriks Transisi Pengukuran H
Matriks
transisi
pengukuran
H
diturunkan
dari
fungsi
pengukuran
h(x ) terhadap variabel keadaan x , dimana persamaan pengukuran yang
dipakai adalah :
⎡ x m ,k ⎤ ⎡ x k ⎤ ⎡ xk ⎤ ⎢ y ⎥ = ⎢ ⎥ Æ h(x ) = ⎢ ⎥ ⎣ yk ⎦ ⎣ m ,k ⎦ ⎣ y k ⎦
(C.5)
maka :
⎡ ∂h1 ∂h(x ) ⎢ ∂x1 =⎢ H= ∂x ⎢ ∂h 2 ⎢⎣ ∂x1
∂h 1 ∂x 2 ∂h 2 ∂x 2
∂h1 ∂x 3 ∂h 2 ∂x 3
∂h1 ⎤ ∂x 4 ⎥ ⎥ ∂h 2 ⎥ ∂x 4 ⎥⎦
(C.6)
Sehingga didapatkan :
⎡1 0 0 0 ⎤ H=⎢ ⎥ ⎣0 0 1 0⎦
(C.7)
Persamaan (C.7) dipakai pada persamaan (3.8) di Bab III.
88