FUSI CITRA SATELIT MULTI-TEMPORAL DENGAN NON-ITERATIF PSEDOPOLAR FOURIER TRANSFORM

FUSI CITRA SATELIT MULTI-TEMPORAL DENGAN NON-ITERATIF PSEDOPOLAR FOURIER TRANSFORM Arya Y. Wijaya1 , Agus Z. Arifin2, Anny Yuniarti3, Wijayanti N. Khotimah4 1,2,3,4

Teknik Informatika, FTIF, Institut Teknologi Sepuluh Nopember Sukolilo, Surabaya 60111 1 [email protected] , 2 [email protected], 3 [email protected], 4 [email protected]

Abstrak Pembuatan peta dari citra satelit biasanya mengalami kesulitan jika hanya memakai sebuah citra. Hal ini dikarenakan seringnya muncul gangguan (misalnya berupa awan) yang berada di atas suatu wilayah. Biasanya, para ahli geografi melakukan fusi citra satelit multi-temporal untuk membuat peta wilayah yang terbebas dari gangguan tersebut. Penelitian ini bertujuan untuk melakukan fusi pada citra satelit multi-temporal untuk mendukung ahli geografi melakukan pemetaan suatu wilayah permukaan bumi secara otomatis. Masukan sistem ini adalah dua citra satelit multi-temporal. Kedua citra satelit tersebut dibawa pada domain frekuensi menggunakan pseudopolar fourier transform untuk dilakukan registrasi. Pada koordinat log polar di domain frekuensi, transformasi geometri citra kedua terhadap citra pertama ditemukan dengan menggunakan phase correlation. Parameter skala perbesaran, sudut rotasi, dan nilai translasi citra kedua terhadap citra pertama yang ditemukan pada tahap registrasi dijadikan dasar pada proses fusi citra. Dengan adanya otomatisasi dalam melakukan fusi citra satelit multi-temporal yang diusulkan dalam penelitian ini, diharapkan pemetaan yang dilakukan lebih objektif, akurat, dan dilakukan dalam waktu relatif singkat. Kata kunci : image fusion, image registration, psudo-polar fourier transform

1.

Pendahuluan

Pembuatan peta dari citra satelit biasanya mengalami kesulitan jika hanya memakai sebuah citra. Hal ini dikarenakan seringnya muncul gangguan (misalnya awan) yang berada di atas suatu wilayah. Biasanya, para ahli geografi memakai lebih dari satu citra satelit pada wilayah yang sama untuk membuat peta wilayah yang terbebas dari gangguan. Citra satelit multi-temporal, yaitu dua citra satelit pada daerah yang sama tetapi diambil pada waktu yang berbeda, digabungkan untuk meningkatkan nilai informasi dari suatu citra satelit. Proses penggabungan yang dilakukan dilakukan secara manual oleh pakar memiliki tingkat subyektivitas yang tinggi sehingga masing-masing pakar bisa memiliki hasil berbeda untuk daerah observasi yang sama. Dengan kemajuan ilmu pada bidang pengolahan citra digital dan visi komputer, metode manual tersebut dapat diganti dengan sebuah metode otomatis. Proses penggabungan, pengkombinasian, penyinergian, atau hal sejenisnya pada dua citra atau lebih menjadi sebuah citra yang memiliki informasi

yang lebih berkualitas dalam ilmu pengolahan citra digital disebut dengan fusi citra [1]. Banyak penelitian tentang fusi yang telah dilakukan. Pada umumnya, penelitian yang dilakukan adalah melakukan fusi citra satelit resolusi tinggi (highresolution panchromatic images/HRPIs) dengan citra resolusi rendah (low-resolution multispectral images/LMRIs). [2,3] menggunakan domain spasial dan melakukan analisis nilai hue dan saturation citra HRPI dan LMRI. (Aiazzi, 2002) menggunakan highpass filtering untuk melakukan fusi. Selain itu, fungsi wavelet juga dapat digunakan untuk melakukan fusi [4]. Akan tetapi, hampir semua metode yang ada lebih berfokus pada fusi antara citra yang diambil dari sensor berbeda, bukan pada citra multi-temporal. Sebelum melakukan fusi, hal yang biasanya dilakukan adalah melakukan registrasi citra. Registrasi citra adalah proses menemukan kembali titik-titik yang bersesuaian antara dua citra.. Citra kedua dianggap merupakan versi transformasi geometri dari citra pertama. Registrasi citra pada domain spasial dilakukan dengan cara mencari nilai rata-rata, median, atau ukuran statistika lainnya pada

Gambar 1. Algoritma Fusi Citra Satelit Multi-Temporal dengan Non-Iteratif Psedopolar Fourier Transform

parameter yang didapatkan pada proses registrasi akan digunakan dasar dalam fusi citra satelit multitemporal. 2. a

b

c

Gambar 2. Pseudopolar-grid (contoh: input citra 8x8) a.P1 b.P2 c. P = P1  P2

setiap nilai derajat keabuan (grayscale) atau RGB citra [5]. Registrasi citra pada domain spasial bekerja dengan baik ketika diaplikasikan terhadap citra yang memiliki tingkat ketidakteraturan kecil. Registrasi citra pada domain frekuensi dilakukan dengan memanfaatkan properti Transformasi Fourier. Metode berbasis Transformasi Fourier mampu memperkirakan skala perbesaran, rotasi, dan translasi lebih akurat dibandingkan dengan metode pada domain spasial. Sebagian besar pendekatan yang berdasarkan Transformasi Fourier memanfaatkan shift property Transformasi Fourier yang menyediakan perkiraan translasi yang akurat dengan phase correlation [6]. Metode yang diusulkan pada penelitian ini adalah melakukan fusi citra satelit multi-temporal dengan terlebih dahulu melakukan registrasi citra berbasis Pseudo-polar Fourier Transform (PPFT) [7]. Akan tetapi, regitrasi citra berbasis PPFT pada [7] menggunakan iterasi agar mendapatkan hasil optimal. Pada penelitian ini, iterasi akan dihilangkan dan diganti dengan metode adaptif untuk mendapatkan estimasi skala perbesaran, sudut rotasi, dan nilai translasi yang optimal. Paramater-

Metodologi Algoritma yang dibangun meliputi tiga tahapan utama sebagaimana digambarkan dalam Gambar 1, yaitu Tahap Representasi PPFT, Tahap Estimasi dan Tahap Fusi. Tahap Representasi PPFT bertujuan untuk membawa citra dari domain spasial menjadi domain frekuensi. Selanjutnya Tahap Estimasi bertujuan untuk mendapatkan nilai transformasi geometri (perbesaran, rotasi, dan pergeseran) citra kedua terhadap citra pertama. Terakhir, Tahap Fusi bertujuan untuk menggabungkan citra kedua yang telah diperbesar, dirotasi dan digeser sesuai hasil registrasi kepada citra pertama. 2.1 Tahap Representasi PPFT Definisi PPFT Pseudo-polar Fourier Transform (PPFT) adalah Transformasi Fourier 2D dari sebuah citra yang letakkan di atas pseudopolar-grid. Secara lengkap, pseudopolar-grid diberikan oleh himpunan (1) P  P1  P2 dimana  2l N   N P1    k , k  |   l  , N  k  N  N 2 2    

(2)

 2l  N N  P2   k , k  |   l  , N  k  N  . N 2 2   

(3)

dan

Untuk mengilustrasikan himpunan P1 dan P2 dapat dilihat pada Gambar 2(a) dan Gambar 2(b). Pseudopolar-grid P diilustrasikan pada gambar 2(c). Pada Gambar 2(a) dan 2(b), k menyajikan nilai yang disebut sebagai pseudoradius dan l menyajikan nilai yang disebut sebagai pseudoangle. Resolusi dari pseudopolar-grid adalah N+1 dalam bagian angular dan M=2N+1 pada bagian radial. Dengan menggunakan representasi (r,θ), pseudopolar-grid diberikan dengan (4) P1 ( k , l )  ( rk1 , l1 ) ; P2 ( k , l )  (rk2 , l2 ), 2

2

 l   l  rk1  k 4   1; rk2  k 4   1 N N

(5)

 2l   2l   l1   / 2  arctan  ;  l2  arctan  N N

(6)

dimana k=-N,....,N dan l = -N/2,....,N/2. PPFT definisikan sebagai sampel dari Transformasi Fourier yang diberikan pada (2.1) di atas pseudopolar-grid P yang diberikan pada (2.10). j Secara detail, PPFT Iˆpp ( j  1, 2) adalah sebuah transformasi linear, dimana terdefinisi untuk k=N,....,N dan l=-N/2,....,N/2, sebagai  2l  Iˆ1pp  Iˆ  k , k   N   2  2i  I (u, v ) exp    ku  kv    M N 2 l    Iˆ k , k  N  

Iˆ1pp 

(7)

N / 2 1



u ,v  N / 2

2 Iˆ pp

(8)

N / 2 1

 2  2i  2 Iˆ pp   I (u, v) exp   ku  kv  N  u ,v   N / 2  M

dimana Iˆ adalah I pada domain frekuensi. Sebagaimana dapat dilihat pada Gambar 2(c), untuk setiap sudut yang telah ditentukan sebesar l, sampel dari pseudopolar-grid memiliki space yang sama pada bagian radial. Akan tetapi, space ini berbeda untuk sudut yang berbeda. Demikian pula, grid memiliki space yang tidak sama dalam bagian angular, tetapi memiliki space kemiringan yang sama. Secara detail, 2 (9)  tan  1pp (l )  cot l11  cot l1  N 2 2  tan pp (l )  cot l21  cot  l2  N

(10)

2 dimana  1pp dan  pp diberikan pada (6).

selanjutnya adalah melakukan penghitungan nilai translasi sepanjang sumbu logr dan sumbu θ. Penghitungan pergeseran ini dilakukan dengan metode phase correlation.

Properti penting PPFT adalah bahwa transformasi ini memiliki kemampuan invert. Selain itu, PPFT forward dan invert dapat diaplikasikan dengan sebuah komputasi yang cepat dengan bantuan FrFT. Dan yang lebih penting lagi, algoritma ini tidak membutuhkan regriding atau interpolasi sehingga memliki keakuratan yang tinggi. Fractional Fourier Fourier (FrFT) Kompleksitas penghitungan PPFT dapat ditekan dengan bantuan FrFT. FrFT adalah algoritma cepat dengan komputasi O(N logN) yang dapat memetakan Transformasi Fourier Diskrit (DFT) di atas beberapa himpunan dari N titik pada sebuah keliling lingkaran. Lebih spesifik, diberikan sebuah vektor C dengan panjang N+1, C = (C(u), u = -N/2,......,N/2),   R . FrFT didefinisikan sebagai : N /2

( FN1C )(k ) 

 C (u) exp[2iku /( N  1)];

u N / 2

k   N / 2,..., N / 2

(11)

Algoritma Tahap Representasi PPFT 1. Input citra I1 dan I2, jika citra memiliki ukuran yang berbeda, maka dilakukan zeropadding (memberikan nilai 0 pada piksel tambahan) sehingga kedua cira memiliki ukuran n x n. 2. Hitung magnitude PPFT pada masing-masing citra I1 dan I2 . 3. Lakukan interpolasi dari Pseudopolar-grid menjadi Log-polar grid sebagaimana digambarkan pada Gambar 3(a). Titik-titik yang berwarna abu-abu adalah elemen pada pseudopolar-grid. Sedangkan titik-titik warna hitam adalah elemen dari pada log-polar grid. Nilai sebuah titik hitam akan mengambil nilai sebuah titik abu-abu yang berada paling dekat dengannya. 4. Dikarenakan nilai interval yang sangat besar (01014) pada setiap nilai yang berada pada koordinat log-polar grid, lakukan operasi log pada seluruh elemen sehingga memiliki interval (0-14). Agar dapat diamati dan memiliki tingkat kontras yang baik, interpolasikan nilai tersebut pada interval grayscale (0-255) sehingga memiliki hasil seperti terlihat pada Gambar 3(b). 2.2 Tahap Estimasi Tahap Estimasi bertujuan untuk mendapatkan nilai transformasi geometri (perbesaran, rotasi, dan pergeseran) citra ke-2 (I2) terhadap citra ke-1 (I1). Setelah mendapatkan nilai PPFT pada log-polar grid. Phase Correlation Apabila I2 adalah citra I1 yang mengalami pergeseran sebesar (Δx, Δy) sehingga, (12) I 1 ( x, y )  I 2 ( x   x, y   y )

Selanjutnya dilakukan penghitungan sebagaimana formula dibawah ini,

Sumbu logr

Sumbu θ

R

a

b

Gambar 3. a. Interpolasi dari Pseudopolar-grid (abuabu) ke polar-grid (hitam). b. Representasi grayscale PPFT yang telah diinterpolasikan pada polar-grid.

a

Iˆ1 Iˆ2 * . | Iˆ1Iˆ2 * |

R (14)

I* adalah complex conjugate dari I. Selanjutnya dicari phase correlation r pada domain spasial dimana r adalah (15) r  1{R} . Nilai (Δx, Δy) dapat ditemukan dengan mencari letak puncak dari dari r yaitu (16) (x, y )  arg max{r} . ( x, y )

Algoritma Tahap Estimasi 1. Dengan menggunakan phase correlation, lakukan penghitungan pergeseran nilai PPFT I2 terhadap I1 yang telah dibawa pada koordinat log-polar pada tahap sebelumnya, sehingga didapatkan nilai pergeseran (Δlogr , Δθ). 2. Nilai perbesaran citra input I2 terhadap I1 sebesar s dapat dicari dari nilai Δlogr [7]. Apabila sumbu logr terdiri dari q buah piksel maka skala perbesaran relatif citra I2 terhadap I1 sebesar s didapatkan dengan

b

 max(log r )    log r  q 1  

3. c

d

Gambar 4. a. Citra “landsat” (I1). b. Citra I2. c. PPFT citra I1 dan I2 pada koordinat log-polar yang telah di bawa pada nilai grayscale. d. Hasil Fusi.

 t 

Tabel 1. Hasil uji coba akurasi registrasi Estimasi Rotasi (derajat)

Skala

Rotasi (derajat)

Error Skala

Error Rotasi (derajat)

1,1

115

1,10

114,96

0,00

0,04

1,5

150

1,49

150,12

0,01

0,12

2,0

246

2,02

246,09

0,02

0,09

3,0

258

3,02

258,05

0,02

0,05

4,0

044

4,01

044,30

0,01

0,30

5,0

003

5,04

003,16

0,04

0,16

Skala

maka besar pergeseran I2 terhadap I1 sebesar (Δx, Δy) dapat ditemukan secara akurat dengan menggunakan phase correlation [8]. Metode estimasi pergeseran dengan menggunakan phase correlation dapat dimulai dengan melakukan Transformasi Fourier 2D pada I1 dan I2 sehingga secara berturut-turut menghasilkan Iˆ1 dan Iˆ2 . (13) Iˆ  {I } ; Iˆ  {I } 1

1

2

2

s  base , (17) dimana base adalah basis logaritma yang digunakan pada saat interpolasi ke log-polar grid pada Tahap Representasi PPFT. Nilai rotasi citra input I2 terhadap I1 sebesar θ dapat dicari dari nilai Δθ [7]. Apabila sumbu θ terdiri dari t buah piksel maka sudut rotasi antara citra I2 terhadap I1 sebesar θ dapat ditemukan dengan  180  . (18)    

4.

5.

Nilai θ yang ditemukan memiliki ambiguitas yaitu θ atau θ+π [7]. Ambiguitas dari θ diselesaikan dengan menggunakan 2D phase correlation dimana citra I1 diputar sebesar θ dan θ+π dan di hitung phase correlationnya terhadap I2. Nilai phase correlation terbesar diantara kedua pasang citra tersebut merupakan sudut θ yang sebenarnya. Optimasi (smoothing) nilai skala dan rotasi yang didapat dilakukan dengan pengambilan q buah nilai s dari tetangga skala estimasi dan q buah nilai θ dari tetangga sudut estimasi. Dilakukan phase correlation dari q+1 buah kandidat nilai skala dan sudut rotasi yang menghasilkan perbesaran dan rotasi paling akurat. Nilai skala perbesaran (s) dan sudut rotasi (θ) yang paling akurat adalah yang menghasilkan nilai phase correlation yang maksimum. Pergeseran/translasi I2 terhadap I1 sebesar (Δx, Δy) ditemukan dengan melakukan phase correlation I1 terhadap citra I2 yang telah di perbesar sebesar 1/s dan diputar sebesar –θ.

a. Indonesia

b. Singapura

c. Malaysia

d. Thailand Gambar 4. Fusi 4 pasang Citra Satelit Multi-temporal. Dari kiri ke kanan: citra input I1 , citra input I2, Hasil rekonstruksi citra I2 sesuai dengan nilai estimasi yang ditemukan pada registrasi, Hasil Fusi Citra Satelit Multi-temporal

2.3 Tahap Fusi Tahap Fusi bertujuan untuk menggabungkan citra kedua yang telah diperbesar, dirotasi dan digeser sesuai hasil registrasi kepada citra pertama. Jika IF adalah citra hasil fusi I1 dengan I2 yang telah diketahui nilai skala perbesaran sebesar s, rotasi sebesar θ, dan translasi sebesar (Δx, Δy) terhadap I1, maka (19) I F  5( I1  I 2 *) / 8 ,

dimana I2* adalah citra I2 yang telah diperbesar dengan 1/s, diputar sejauh –θ dan digeser sejauh (Δx, -Δy) 3.

Uji Coba

3.1 Uji Coba Akurasi Uji coba ini bertujuan untuk menentukan akurasi estimasi skala perbesaran dan sudut rotasi oleh algoritma yang diusulkan. Citra uji I1 adalah citra “landsat” dan citra I2 adalah citra “landsat”

Tabel 2. Hasil uji coba pada citra multi-temporal No

Nama Pasangan Citra

Transformasi Geometri I2 terhadap I1 Skala

Sudut Rotasi

Translasi

1

Indonesia

1.10

2.02

(6, -5)

2

Singapura

1.00

0.00

(35,-1)

3

Malaysia

1.20

4.50

4

Thailand

1.05

-174.88

(8,21) (15,-17)

yang telah dilakukan perbesaran sebesar s dan rotasi sebesar θ. Salah satu contoh uji coba dapat dilihat pada Gambar 4. Secara berutan: Gambar 4(a) adalah citra inputI1, Gambar 4(b) adalah citra input I2, Gambar 4(c) adalah PPFT citra I1 dan I2 pada koordinat logpolar yang telah di bawa pada nilai grayscale, dan Gambar 4(d) hasil fusi citra. Secara keseluruhan, hasil uji coba dapat dilihat pada Tabel 1. Terlihat bahwa algoritma yang digunakan memiliki tingkat kesalahan yang rendah dalam menemukan skala dan sudut estimasi. 3.2 Uji Coba pada Citra Satelit Multi-temporal Uji coba ini bertujuan untuk melihat seberapa jauh algoritma yang diusulkan mampu melakukan fusi dua citra satelit multi-temporal. Data uji citra satelit pada penelitian ini adalah terdiri dari empat pasang citra satelit multi-temporal yang didapatkan dari Aplikasi Google Earth. Citra multi-temporal adalah citra dengan objek yang sama tetapi diambil pada waktu yang berlainan. Citra pertama I1 adalah citra yang berawan (gangguan) dan lebih baru tanggal pengambilannya dibanding citra kedua I2. Dimensi seluruh citra tersebut masing-masing adalah 580x580. Hasil uji coba empat pasang citra satelit multitemporal secara visual dapat dilihat pada Gambar 5. Secara berurutan, dari kiri ke kanan adalah: citra I1, citra I2, citra I2 yang di diperbesar - diputar - digeser sesuai hasil regitrasi, dan paling kanan adalah citra hasil Fusi. Secara kuantitatif, hasil estimasi skala perbesaran, sudut rotasi dan pergeseran dapat dilihat pada Tabel 2. Terlihat bahwa tingkat kesalahan (error) tidak dapat dihitung secara kuantitatif karena groundthruth tidak diketahui. Berbeda dengan uji coba 3.1 yang nilai perbesaran dan rotasi diketahui sebelumnya, sehingga tingkat kesalahan dapat ditemukan. Untuk menyatkan hasil fusi pada uji coba ini hanya dapat dinyatakan secara kuantitatif bagus, atau tidak bagus dari hasil pengamatan. Dari pengamatan visual dapat dinyatakan bahwa algoritma yang diusulkan memiliki hasil fusi yang bagus. Hal ini dikarenakan pada citra fusi terlihat bahwa objek-objek terlihat tetap terhubung walau

citra fusi merupakan gabungan dari dua citra satelit multi-temporal. 4.

Kesimpulan Dari hasil ujicoba yang dilakukan dapat disimpulkan bahwa algoritma yang diusulkan memiliki akurasi yang tinggi dalam melakukan estimasi perbesaran, rotasi dan pergeseran. Ketika algoritma yang diusulkan diterapkan pada citra satelit multi-temporer, dapat dilihat bahwa hasil fusi yang dilakukan memiliki hasil yang bagus. Kedepannya, perlu dikembangkan metode segmentasi hambatan yang terdapat pada citra satelit pertama (misalnya awan). Sebelum dilakukan fusi, hambatan tersebut dihapus sehingga bagian detail yang hilang/dihapus tersebut dapat diambilkan dari citra satelit yang kedua pada koordinat yang bersesuian. 5. Daftar Pustaka [1] Zhijun Wang, dkk., 2005, “A Comparative Analysis of Image Fusion Methods,” IEEE Trans. Geosci. Remote Sens.,vol. 43, no. 6, pp. 81–84, Jun. 2005. [2] Liu,J.G., 2000 , “Smoothing filter-based intensity modulation: A spectral preserve image fusion technique for improving spatial details,” Int. J. Remote Sens., vol. 21, no. 18, pp. 3461–3472. [3] Tu, T.M. dkk.,2001, “A new look at IHS-like image fusion methods,” Inf. Fusion, vol. 2, no. 3, pp. 177–186. [4] Ranchin, T. , 2003, “Image fusion—The ARSIS concept and some successful implementation schemes,” ISPRS J. hotogramm. Remote Sens., vol. 58, pp. 4–18, 2003. [5] Wolberg, G. dan Zokai. (2000), “Robust image registration using log-polar transform”, Proc. IEEE Int. Conf. Image Processing, Vancouver, BC, Canada [6] Reddy, S. dan Chatterji, B.N, 1996, “An FFTBased Technique for Translation, Rotation, and Scale-Invariant Image Registration”, IEEE Trans. Image Processing, Vol. 3, No. 8, hal. 1266-1270. [7] Keller, Yosi dkk., 2005, “Pseudopolar-Based Estimation of Large Translations, Rotation”, IEEE Transactions on Image Processing, Vol. 14, No.1, Hal. 12-22. [8] Kuglin, C.D. and Hines, D.C., 1975, “The phase correlation image alignment method,” in Proc. IEEE Conf. Cybernetics and Society, Hal. 163–165.