Pengolahan Citra di Kawasan Frekuensi

BAB 6 Pengolahan Citra di Kawasan Frekuensi

Setelah bab ini berakhir, diharapkan pembaca mendapatkan berbagai pengetahuan berikut dan mampu mempraktikkannya.  Pengolahan citra di kawasan spasial dan kawasan frekuensi  Alihragam Fourier  Fourier 1-D  Fourier 2-D  Fast Fourier transform  Visualisasi pemrosesan FFT  Penapisan pada kawasan frekuensi  Filter lolos-rendah  Filter lolos-tinggi  Penapisan dengan pendekatan high frequency emphasis

168

Pengolahan Citra, Teori dan Aplikasi

6.1 Pengolahan Citra di Kawasan Spasial dan Kawasan Frekuensi Citra dapat ditransformasikan di kawasan spasial maupun di kawasan frekuensi. Dua cara untuk melakukan transformasi citra ditunjukkan pada Gambar 6.2. Sejumlah contoh transformasi di kawasan spasial atau keruangan telah dibahas di Bab 5. Namun, sebagai alternatif citra perlu diproses di kawasan frekuensi agar penentuan daerah frekuensinya dapat lebih ketat dan tepat. Untuk itu diperlukan pasangan transformasi dan transformasi-balik sebelum dan sesudah penapisan.

Gambar 6.1 Transformasi citra dapat diproses melalui kawasan spasial maupun frekuensi



Dalam bahasa Indonesia, istilah lain yang identik dengan transformasi adalah alihragam.



Adanya pasangan alihragam dan alihragam-balik tentu saja menambah beban komputasi.

Salah satu alihragam yang biasa dipakai di kawasan frekuensi adalah alihragam

Fourier.

bermunculan

pula

Sejak

algoritma

alihragam

Fourier

macam-macam alihragam yang lain,

ditemukan,

telah

seperti transformasi

Pengolahan Citra di Kawasan Frekuensi

169

gelombang-singkat (wavelet), transformasi Radon,

dan DCT (Discrete Cosine

Transform). 6.2 Alihragam Fourier Alihragam Fourier (Fourier transform) merupakan salah satu jenis alihragam ke kawasan frekuensi yang banyak dipakai pada pengolahan citra. Alihragam ini dimanfaatkan untuk memetakan citra dari kawasan spasial ke dalam kawasan frekuensi. Disamping untuk melihat karakteristik spektrum citra, juga menjadi

bagian

pemrosesannya.

Citra

dapat

diamati

sebagai

kumpulan

gelombang sinusoid dengan frekuensi, amplitudo, dan fase yang berbeda-beda. Meskipun pada zaman sekarang terdapat berbagai alihragam sebagai alternatif alihragam

Fourier, konsep yang mendasari alihragam Fourier perlu dimengerti.

Lagipula, beberapa pemrosesan masih bertumpu pada alihragam Fourier. Berdasarkan temuan ahli fisika dari Prancis bernama Baptiste Joseph Fourier (1768-1830), semua fungsi yang bersifat periodis, betapapun kompleks fungsi tersebut,

dapat

dinyatakan sebagai penjumlahan sinusoid.

Kuncinya

terletak pada komposisi amplitude dan fase sinus setiap frekuensi. Begitu pula pada citra. Sebagai gambaran, suatu isyarat berdimensi satu pada Gambar 6.1(a) dapat disusun atas tiga gelombang sinusoid seperti terlihat pada Gambar 6.1(b).

Gambar 6.1 Contoh isyarat yang tersusun atas tiga sinusoid

170


6.3 Fourier 1-D Penerapan Discrete Fourier Transform (DFT) atau alihragam Fourier diskret pada citra berdimensi satu disajikan pada pembahasan berikut. Misalnya, terdapat fungsi f(x) yang terdiri atas N data (f(0), f(1), f(2), f(3), f(4), …, f(N-1)). Jika dikenakan DFT, akan diperoleh hasil alihragam berupa F(u) berupa F(u) = (F(0), F(1), F(2), F(3), F(4), .., F(N-1)) Perhatikan bahwa jumlah data diskret N yang sama di kawasan frekuensi, yang sejalan dengan hukum kelestarian informasi. F(u) diperoleh melalui persamaan: 𝐹 (𝑢) =

1 𝑁

2𝜋𝑢𝑥 2𝜋𝑢𝑥 ∑𝑁−1 ] − 𝑗 𝑠𝑖𝑛 [ ]) 𝑥=0 𝑓(𝑥)(𝑐𝑜𝑠 [ 𝑁

𝑁

(6.1)

atau 1

𝐹 (𝑢) = 𝑁 ∑𝑁−1 𝑥=0 𝑓( 𝑥 ) 𝑒𝑥𝑝 [−𝑗

2𝜋𝑢𝑥 𝑁

] , dengan u = 0,1,2,…,N-1

(6.2)

Pada rumus di depan, j menyatakan √ −1. Dengan demikian, hasil transformasi Fourier berupa bilangan kompleks. Adapun alihragam-baliknya berupa: 2𝜋𝑢𝑥 𝑓 (𝑥 ) = ∑𝑁−1 𝑥=0 𝐹 ( 𝑢) 𝑒𝑥𝑝[𝑗 𝑁 ] , dengan u = 0,1,2,…,N-1

(6.3)

Sebagai contoh, terdapat f(x) = (2, 4, 1, 5). Alihragam Fourier-nya seperti berikut. 1

𝐹 (0) = 4 ∑3𝑥=0 𝑓(𝑥)(𝑐𝑜𝑠 [

2𝜋0𝑥 4

] − 𝑗 𝑠𝑖𝑛 [

2𝜋0𝑥 4

])

= (f(0)(cos(0)-j sin(0)) + f(1)(cos(0)-j sin(0)) + f(2)(cos(0)-j sin(0)) + f(3)( (cos(0)-j sin(0))) / 4 = (f(0) + f(1) + f(2) + f(3)) / 4 = (2 + 4 + 1 + 5) / 4 = 12 / 4 = 3 1 2𝜋1𝑥 2𝜋1𝑥 𝐹 (1) = 4 ∑3𝑥=0 𝑓(𝑥)(𝑐𝑜𝑠 [ 4 ] − 𝑗 𝑠𝑖𝑛 [ 4 ])


171

= (f(0)(cos(0)-j sin(0)) + f(1)(cos(𝜋/2)-j sin(𝜋/2)) + f(2)(cos(𝜋)-j sin(𝜋)) + f(3)( (cos(3𝜋/2)-j sin(3𝜋/2))) / 4 = (2 (1-0) + 4(0-j) + 1(-1-0) + 5(0+j)) / 4 = (1+j)/4 = 0,25 + j0,25 1

𝐹 (2) = 4 ∑3𝑥=0 𝑓(𝑥)(𝑐𝑜𝑠 [

2𝜋2𝑥 4

] − 𝑗 𝑠𝑖𝑛 [

2𝜋2𝑥 4

])

= (f(0)(cos(0)-j sin(0)) + f(1)(cos(𝜋)-j sin(𝜋)) + f(2)(cos(2𝜋)-j sin(2𝜋)) + f(3)( (cos(3𝜋)-j sin(3𝜋))) / 4 = (2 (1-0) + 4(-1-0) + 1(1-0) + 5(-1-0) ) / 4 = -6 / 4 = -1,50

1 2𝜋3𝑥 2𝜋3𝑥 𝐹 (3) = 4 ∑3𝑥=0 𝑓(𝑥)(𝑐𝑜𝑠 [ 4 ] − 𝑗 𝑠𝑖𝑛 [ 4 ])

= (f(0)(cos(0)-j sin(0)) + f(1)(cos(3𝜋/2)-j sin(3𝜋/2)) + f(2)(cos(3𝜋)-j sin(3𝜋)) + f(3)( (cos(9𝜋/2)-j sin(9𝜋/2))) / 4 = (2 (1-0) + 4(0+j) + 1(-1-0) + 5(0-j) ) / 4 = (1 - j)/4 = 0,25 - j0,25

Gambar 6.2 memperlihatkan citra asli dan hasil transformasi Fourier.

172


Gambar 6.2 Citra dimensi satu dan hasil transformasi Fourier Pada Gambar 6.2, DFT-1 menyatakan transformasi balik dari kawasan frekuensi ke kawasan spasial. Perhatikan bahwa data asli f(x) hanya 4, tetapi hasil alihragam ada 8, seolah ada tambahan informasi. Sekarang akan ditunjukkan pelaksanaan alihragam-baliknya. Perhatikan cara menghitungnya. 𝑓 (0) = ∑3𝑥=0 𝐹(𝑥)(𝑐𝑜𝑠 [

2𝜋0𝑥 4

] + 𝑗 𝑠𝑖𝑛 [

2𝜋0𝑥 4

])

= F(0)(cos(0)+j sin(0)) + F(1)(cos(0)+j sin(0)) + F(2)(cos(0)+j sin(0)) + F(3)( (cos(0)+j sin(0))) / 4 = F(0) + F(1) + F(2) + F(3) = 3+ j0 + 0,25 + j0,25 - 1,5+j0 + 0,25-j0,25 =2 3

2𝜋1𝑥 2𝜋1𝑥 𝑓(1) = ∑ 𝐹(𝑥)(𝑐𝑜𝑠 [ ] + 𝑗 𝑠𝑖𝑛 [ ]) 4 4 𝑥=0

= F(0)(cos(0)+j sin(0)) + F(1)(cos(𝜋/2)+j sin(𝜋/2)) + F(2)(cos(𝜋)+j sin(𝜋)) + F(3)( (cos(3𝜋/2)+j sin(3𝜋/2)) = 3 (1+0) + (0,25 + j0,25)(0+j) - 1,5 (-1+0) + (0,25-j0,25)(0-j)


173

= 3 + j0,25 -0,25 + 1,5 - j0,25 -0,25 =4 3

2𝜋2𝑥 2𝜋2𝑥 𝑓 (2) = ∑ 𝐹(𝑥)(𝑐𝑜𝑠 [ ] + 𝑗 𝑠𝑖𝑛 [ ]) 4 4 𝑥=0

= F(0)(cos(0)-j sin(0)) + F(1)(cos(𝜋)-j sin(𝜋)) + F(2)(cos(2𝜋)-j sin(2𝜋)) + F(3)(cos(3𝜋)-j sin(3𝜋)) = 3 (1-0) + (0,25 + j0,25) (-1-0) - 1,5 (1-0) + ((0,25-j0,25))(-1-0) = 3 – 0,25 –j0,25 -1,5 – 0,25 +j0,25 =1 3

2𝜋3𝑥 2𝜋3𝑥 𝑓(3) = ∑ 𝐹(𝑥)(𝑐𝑜𝑠 [ ] + 𝑗 𝑠𝑖𝑛 [ ]) 4 4 𝑥=0

= (F(0)(cos(0)+j sin(0)) + F(1)(cos(3𝜋/2)+j sin(3𝜋/2)) + F(2)(cos(3𝜋)+j sin(3𝜋)) + F(3)( (cos(9𝜋/2)+j sin(9𝜋/2))) / 4 = 3 (1-0) + (0,25 + j0,25) (0-j) - 1,5 (-1-0) + (0,25-j0,25)(0+j) = 3 - j0,25 + 0,25 + 1,5 + j0,25 + 0,25 = 5

Tampak bahwa f(x) tidak mengandung komponen imajiner seperti aslinya. Berikut adalah contoh fungsi yang digunakan untuk menghitung DFT berdimensi satu. Program : dft1d.m

174


function [Re, Im] = dft1d(Fx) % DFT1D Digunakan untuk memperoleh DFT dimensi satu. % Hasil: Re berisi bagian real dan % Im berisi bagian imajiner n = length(Fx); % Jumlah nilai dalam fungsi Fx for u = 0 : n - 1 Re(u+1) = 0; Im(u+1) = 0; for x = 0 : n - 1 radian = 2 * pi * u * x / n; cosr = cos(radian); sinr = -sin(radian); Re(u+1) = Re(u+1) + Fx(x+1) * cosr; Im(u+1) = Im(u+1) + Fx(x+1) * sinr; end Re(u+1) = Re(u+1) / n; Im(u+1) = Im(u+1) / n; end

Akhir Program

Contoh penggunaan fungsi dft1d

ditunjukkan di bawah ini. Perhatikan Fx

perlu ditranspos. >> Fx = [2,4,1,5]';  >> [Re,Im] = dft1d(Fx)  Re = 3.00000

0.25000

-1.50000

0.25000

0.25000

-0.00000

-0.25000

Im = 0.00000 >> Hasil di atas sesuai dengan perhitungan manual di depan. Alihragam-baliknya

dapat

berikut.

Program : idft1d.m

diperoleh

melalui fungsi idft1d

seperti


175

function Fx = idft1d(Fu) % IDFT1D Digunakan untuk melaksanakan transformasi balik % 1D DFT. % Masukan: Fu berupa bilangan kompleks n = length(Fu); % Jumlah nilai dalam fungsi Fu for u = 0 : n - 1 Fx(u+1) = 0; for x = 0 : n - 1 radian = 2 * pi * u * x / n; cosr = cos(radian); sinr = sin(radian); Fx(u+1) = Fx(u+1) + Fu(x+1) * (cosr+ j*sinr); end end Fx = real(Fx); % Peroleh bagian real

Akhir Program

Berikut adalah contoh penggunaan DCT dan alihragam-baliknya: >> Fx = [2,4,1,5]';  >> [Re,Im] = dft1d(Fx);  >> F = idft1d(Re+Im*j)  F = 2.0000

4.0000

1.0000

5.0000

>> Perhatikan, argumen pada fungsi idft1d

berupa bilangan kompleks, tetapi

tidak ada komponen imjinernya atau nilai imajinernya nol. Pengalian dengan j pada Im*j dengan isi Fx.

dipakai untuk membentuk bagian imajiner. Hasilnya terlihat sama

176


6.4 Fourier 2-D Suatu citra diskret berdimensi dua f(x, y) dapat dinyatakan sebagai deret Fourier, yang dituliskan seperti berikut: 𝑢𝑥

−1 𝑁−1 𝐹(𝑢, 𝑣) = ∑𝑀 𝑦=0 ∑𝑥 =0 𝑓(𝑥, 𝑦)(cos(2𝜋 ( 𝑁 +

𝑣𝑦 𝑀

𝑢𝑥

)) − 𝑗 sin (2𝜋 ( 𝑁 +

𝑣𝑦 𝑀

))) (6.4)

Dalam hal ini, citra berukuran MxN (M baris dan N kolom). Komponen v bernilai dari 0 sampai dengan M-1 dan u bernilai dari 0 sampai dengan N-1. Dalam hal ini, u dan v menyatakan frekuensi, sedangkan nilai F(u, v) dinamakan koefisien Fourier atau spektrum frekuensi diskret. Adapun alihragam-baliknya berupa: 𝑓(𝑦, 𝑥) =

1 𝑀𝑁

𝑢𝑥

𝑁−1 ∑𝑀−1 𝑦=0 ∑ 𝑥=0 𝐹(𝑣, 𝑢)(cos (2𝜋 (

𝑁

+

𝑣𝑦 𝑀

𝑢𝑥

)) + 𝑗 sin (2𝜋 (

𝑁

+

𝑣𝑦 𝑀

)))

(6.5)

Berdasarkan rumus di atas, setiap piksel akan ditransformasikan dengan kompleksitas berupa O(MN) atau O(N 2 ) jika ukuran citra berupa NxN. Dengan demikian, untuk citra berukuran NxN, kompleksitas berupa O(N 4 ). Tentu saja, untuk ukuran yang besar akan diperlukan waktu yang sangat lama. Itulah sebabnya, dalam praktik cara tersebut dihindari karena terdapat metode lain yang jauh lebih cepat dalam melakukan transformasi. Sekedar untuk kepentingan ilustrasi, alihragam Fourier yang menggunakan Persamaan 6.4 dapat diimplementasikan seperti berikut. Program : dft2d.m

function [Re, Im] = dft2d(berkas) % DFT2D Digunakan untuk memperoleh DFT dimensi dua. % Masukan: nama berkas berskala keabuan % Hasil: Re berisi bagian real dan % Im berisi bagian imajiner Fx = double(imread(berkas)); [m, n] = size(Fx); % Ukuran citra % m = jumlah baris % n = jumlah kolom for v = 0 : m -1 for u = 0 : n - 1


Re(v+1, Im(v+1, for y = for

u+1) = 0; u+1) = 0; 0 : m - 1 x = 0 : n - 1 radian = 2 * pi * … (u * x / n + v * y cosr = cos(radian); sinr = -sin(radian); Re(v+1, u+1) = Re(v+1, u+1) Fx(y+1, x+1) Im(v+1, u+1) = Im(v+1, u+1) Fx(y+1, x+1)

177

/ m); + * + *

... cosr; ... sinr;

end end end end

Akhir Program

Fungsi dft2d digunakan untuk mentransformasikan citra berskala keabuan. Contoh penggunaannya misalnya seperti berikut: >> [Dr, Di] = dft2d(’C:\Image\lena64.png’);  Dengan cara seperti itu, Dr mencatat bagian real dan Di mencatat bagian imajiner. Tabel berikut memberikan contoh waktu pengeksekusian yang dilakukan pada tiga citra dengan ukuran yang berlainan dengan menggunakan fungsi dft2d. Terlihat bahwa semakin besar ukuran citra yang diproses, terjadi peningkatan waktu komputasi yang sangat signifikan. Tabel 6.1 Waktu komputasi menggunakan dft2d Citra

Waktu (detik)

lena64.png (ukuran 64 x 64)

3,51


54,42


870,09

178


Pendekatan lain memungkinkan penghitungan baris dan kolom tidak dilakukan sekaligus. Pertama, transformasi berdimensi satu hanya dikenakan pada baris-baris saja. Selanjutnya, hasilnya ditransformasi menurut kolom. Contoh program dalam bahasa C terdapat di Vandevenne (2007). 6.5 Fast Fourier Transform Suatu metode bernama FFT (Fast Fourier Transform) dibuat untuk mempercepat komputasi alihragam Fourier.

Jika kompleksitas DFT untuk

mentransformasikan sebuah piksel seperti yang tertuang dalam implementasi di depan sebesar O(N 2 ), FFT memiliki kompleksitas sebesar O(N log2 N). Sebagai pembanding, jika N sama dengan 256 maka N 2 sama dengan 65.536, sedangkan N log2 N menghasilkan 256 x 8 atau 2048. Jadi, FFT lebih cepat 32 kali dibandingkan DFT untuk ukuran citra seperti itu. Pada alihragam berdimensi dua, penghematan waktu akan lebih terasa. Peluang adanya penghematan waktu komputasi dapat dilukiskan untuk citra N=8. Nilai 𝑒𝑥𝑝 [𝑗 real cos 𝑒𝑥𝑝 [𝑗

2𝜋𝑢𝑥

2𝜋𝑢𝑥 8

8

] untuk u dan x sama dengan 0,1,2,3,4,5,6,7 dengan nilai

] + 𝑗𝑠𝑖𝑛 𝑒𝑥𝑝 [𝑗

2𝜋𝑢𝑥 8

] ditunjukkan pada Tabel 6.2.

Tabel 6.2 Nilai-nilai 𝑒𝑥𝑝 [𝑗 x

2𝜋𝑢𝑥 8

]

0

1

2

3

4

5

6

7

0

1+j0

1+j0

1+j0

1+j0

1+j0

1+j0

1+j0

1+j0

1

1+j0

1 1 +𝑗 √2 √2

0+j

2

1+j0

0+j

-1+j0

3

1+j0

4

1+j0

1+j0

1+j0

5

1+j0

1 1 +𝑗 √2 √2

0+j

6

1+j0

0+ j

-1+j0

0-j

1+j0

7

1+j0

0+i

1 1 −𝑗 √2 √2

0-j

u

−

−

1 1 −𝑗 2 √ √2

1 1 −𝑗 √2 √2

0+j

−

1 1 +𝑗 √2 √2

0-j +

1 1 −𝑗 2 √ √2

1+j0 −

1 1 +𝑗 √2 √2

-1+j0

−

1+j0 0-j

0+j +

1+j0 -1+j0

1 1 −𝑗 √2 √2

1 1 +𝑗 2 √ √2

1+j0 −

1 1 −𝑗 √2 √2

0+j +

1 1 +𝑗 √2 √2

0-j

−

-1+j0 0-j

0-j −

1+j0 0-j

1 1 +𝑗 2 √ √2

1+j0 −

-1+j0 0-j

1 1 −𝑗 √2 √2

1 1 −𝑗 √2 √2

0-j −

1 1 +𝑗 √2 √2


Terlihat bahwa nilai-nilai 𝑒𝑥𝑝 [𝑗

2𝜋𝑢𝑥 8

179

] banyak yang sama, yaitu 0,

1 √2

, dan 1.

Oleh karena itu, algoritma cepat dibuat untuk tidak mengulang proses perkalian dengan angka-angka yang sama. Bahkan, perkalian dengan angka 1 dapat diloncati karena tidak perlu dilakukan. Demikian pula, perkalian dengan nol sama dengan melompati datanya untuk tidak perlu disertakan dalam perhitungan. Namun, perlu dictata bahwa penghematan proses komputasi ini hanya dapat terjadi untuk jumlah pikel N = 2 n , yaitu 2, 4, 8, 16, dan seterusnya. Untuk citra dengan besar N (dan M) kurang dari angka-angka tersebut dapat dilengkapi dengan piksel-piksel kosong (bernilai intensitas nol). Cara melakukan komputasi dengan FFT dijabarkan oleh Cooley, dkk. (1969). Implementasi dengan Octave berupa fungsi fft dan fft2. Contoh penggunaan fft2 seperti berikut: >> Img = imread(’C:\Image\lena256.png’);  >> F = fft2(Img); 

Pada contoh di atas, citra bernama lena256.png ditransformasikan menggunakan FFT. 6.6 Visualisasi Pemrosesan FFT Sebagaimana telah dibahas di depan, alihragam Fourier menghasilkan bilangan kompleks. Terkait dengan hal itu, terdapat definisi spektrum Fourier seperti berikut: |𝐹 (𝑣, 𝑢)| = √𝑅 2 (𝑣, 𝑢) + 𝐼 2 (𝑣, 𝑢)

(6.6)

dengan R(u,v) menyatakan bagian real dan I(u,v) menyatakan bagian imajiner. Adapun sudut fase transformasi didefinisikan sebagai: 𝐼(𝑣,𝑢) ∅(𝑣, 𝑢) = 𝑡𝑎𝑛 −1 [ ] 𝑅(𝑣,𝑢)

(6.7)

180


Selain itu, terdapat pula istilah power spectrum atau spektrum daya, yang didefinisikan sebagai kuadrat besaran: 𝑃 (𝑣, 𝑢) = |𝐹(𝑣, 𝑢)|2 = 𝑅 2 (𝑣, 𝑢) + 𝐼 2 (𝑣, 𝑢)

(6.8)

Untuk kepentingan analisis secara visual, spektrum dan sudut fase Fourier dapat disajikan dalam bentuk gambar. Berikut adalah contoh untuk melakukan transformasi citra lena256.png dan kemudian menyajikan spektrum. >> Img = imread(’C:\Image\lena256.png’);  >> F = fft2(Img);  >> imshow(abs(F),[]); 

Citra lena256.png dan spektrum Fourier-nya diperlihatkan pada Gambar 6.3(a) dan Gambar 6.5(b).


(a) Citra lena256.png

(c) Visualisasi hasil FFT dengan skala logaritmik intensitas

181

(b) Visualisasi hasil FFT. Hanya di pojok kiri atas yang terlihat putih

(d) Visualisasi hasil FFT. Dengan frekuensi nol ditengahkan

Gambar 6.3 Transformasi Fourier pada citra dan visualisasi

Mengingat nilai dalam spektrum terlalu lebar, penerapan logaritma biasa digunakan hanya untuk kepentingan visualisasi. Sebagai contoh, Gambar 6.3(c) diperoleh melalui: >> S2 = log(1 + abs(F));  >> imshow(S2, []); 

Penambahan angka 1 dimaksudkan untuk menghindari terjadinya log(0).

182


Hasil pada Gambar 6.3(c) menunjukkan keadaan yang seperti berulang yang muncul pada setiap pojok dalam kotak frekuensi. Hal ini disebabkan adanya sifat pengulangan pada transformasi Fourier. Dalam hal ini, nilai pada M/2 menuju ke M-1 adalah pengulangan dari titik asal 0 hingga M/2 – 1. Hal ini juga berlaku pada arah mendatar. Berdasarkan sifat ini,

untuk kepentingan visualisasi, titik

awal (0,0) seringkali diubah agar terletak di tengah-tengah kotak frekuensi, sebagaimana ditunjukkan pada Gambar 6.4.

Gambar 6.4 Frekuensi 0 diletakkan di tengah-tengah kotak frekuensi

Hal seperti itu dapat dikerjakan dengan menggunakan fungsi fftshift yang disediakan oleh Octave. Kegunaan fungsi fftshift adalah untuk mengatur agar komponen

frekuensi

nol

diletakkan

di

tengah-tengah

memberikan gambaran fungsi ini, perhatikan contoh berikut.

X=

spektrum.

Untuk


1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

183

16

Jika dikenakan perintah seperti berikut

fftshift(X)

diperoleh hasil berupa:

ans =

11

12

9

10

15

16

13

3

4

1

2

7

8

5

6

14

Gambar 6.5 menunjukkan model penukaran keempat kuadran yang bersebrangan secara diagonal melalui fftshift.

Gambar 6.5 Penukaran kuadran melalui fftshift

184


Contoh berikut menunjukkan efek penukaran fftshift pada hasil transformasi Fourier: >> G = fftshift(F);  >> imshow(log(1+abs(G)),[])  >>

Hasilnya dapat dilihat pada Gambar 6.3(d). 6.7 Penapisan pada Kawasan Frekuensi Sebagaimana telah diutarakan pada Gambar 6.1, penapisan dapat dilakukan pada kawasan frekuensi. Menurut teorema konvolusi, konvolusi pada kawasan frekuensi dapat dilakukan dengan mengalikan F(v, u) dengan H(v,u) (Gonzalez, dkk., 2004). Dalam hal ini, H(v,u) dinamakan sebagai fungsi transfer filter dan diperoleh melalui pengenaan DFT terhadap h(y,x), yang merupakan kernel konvolusi pada kawasan spasial. Satu hal yang perlu diperhatikan pada kawasan frekuensi, penapisan dapat menimbulkan problem akibat konvolusi. Problem yang dimaksud dikenal dengan nama wraparound error atau spatial aliasing error (Bovik, 2009). Hal ini disebabkan pada kawasan frekuensi terdapat fungsi periodis (yang berulang setelah jarak

tertentu) yang membuat gambar akan diulang (seperti efek

pengubinan) dan akibatnya membuat interferensi pada konvolusi. Gambar 6.6(a) menunjukkan dua citra yaitu f dan h. Adapun Gambar 6.6(b) menunjukkan operasi konvolusi pada koordinat (m, n). Berdasarkan Gambar 6.6(b), piksel (m, n) dihitung sebagai penjumlahan atas perkalian antara piksel di f dan h. Hasilnya tentu saja mengandung kesalahan yaitu ketika perkalian antara f dan h yang berulang ikut dijumlahkan.


185

Gambar 6.6 Problem wraparound error Untuk mengatasi wraparound error, fungsi f dan h dimodifikasi dengan cara memperbesar ukurannya dan bagian yang diperluas diisi dengan nol. Cara perluasan ukuran dan pemberian nilai nol seperti itu dinamakan sebagai zero padding. Gambar 6.7(a) menunjukkan keadaan setelah zero padding dikenakan pada citra f dan h. Adapun Gambar 6.7(b) menunjukkan bahwa konvolusi pada suatu koordinat piksel tidak lagi mengandung wraparound error.

186


Gambar 6.7 Efek zero padding pada konvolusi pengulangan Gonzales, dkk. (2004) menjelaskan cara menentukan ukuran baru untuk citra f dan h. Misalkan, f berukuran a x b dan h berukuran c x d. Kedua fungsi tersebut diperluas menjadi berukuran p x q. Maka, ukuran untuk p dan q harus memenuhi p > a + c -1 dan

(6.9)


187

q > b + d -1

(6.10)

Apabila menggunakan FFT, nilai p dan q dipilih sehingga memenuhi 2r dengan r adalah bilangan bulat. Dalam praktik, nilai p dan q acapkali dibuat sama, yaitu sebesar r. Dengan kata lain, nilai r harus memenuhi persamaan berikut: 2r > a + c -1

(6.11)

2r > b + d -1

(6.12)

dan

Dengan asumsi bahwa ukuran citra f jauh lebih besar daripada citra h, satu pendekatan yang biasa dilakukan untuk memperoleh nilai p dan q untuk kepentingan zero padding menggunakan algoritma seperti berikut. ALGORITMA 6.1 – Menentukan ukuran zero padding Masukan: 

f dan h (masing-masing berupa citra berukuran a x b dan c x d)

Keluaran: 

pxq (ukuran baru citra untuk f dan h)

1. m  max(a, b) 2. Cari nilai r sehingga memenuhi 2r > 2m – 1 3. p 2q 4. q 2q

188


Di Octave dan MATLAB, nilai q berdasarkan Algoritma 6.1 dapat diperoleh melalui: nextpow2(2 * max(a, b))

Adapun contoh berikut menunjukkan pemakaian zero padding. >> Fs = imread('C:\Image\kotatua.png');  >> Hs = fspecial('sobel');  >> size(Fs)  ans = 747

500

>> Ff = fft2(Fs, 512, 512);  >> Hf = fft2(Hs, 512, 512);  >> G = Hf .* Ff;  >> F = real(ifft2(G));  >> Fx = F(1:250, 1:250);  >> imshow(uint8(Fx)) 

Pada contoh di atas, pengaturan 512 didasarkan pada perhitungan untuk zero padding seperti yang telah dibahas di depan, mengingat ukuran citra berupa 250x250. Pada perintah di atas, Fx = F(1:250, 1:250); digunakan untuk mengambil citra seperti ukuran citra semula. Gambar 6.8 menunjukkan contoh citra yang diproses dan hasil pemfilteran.


189

Gambar 6.8 Contoh pemfilteran dengan filter Sobel pada kawasan frekuensi

Pada perintah di depan Hs = fspecial('sobel'); digunakan untuk memperoleh filter Sobel. Fungsi fspecial disediakan oleh paket Image Processing. Implementasi fungsi Sobel yang tidak menggunakan fspecial akan dibahas pada Bab 10.

Skrip berikut ditujukan untuk memudahkan dalam menangani penapisan pada kawasan frekuensi melalui FFT. Program : filterdft.m

function F = filterdft(berkas, H) % FILTERDFT Digunakan untuk melaksanakan pemfilteran % pada kawasan frekuensi menggunakan FFT. % Masukan: % berkas - nama citra % H - kernel pada kawasan spasial % Keluaran: % F - citra yang telah difilter

190


Fs = double(imread(berkas)); [a, b] = size(Fs); %Peroleh ukuran citra % r p q

Menentukan ukuran baru untuk perluasan citra = nextpow2(2 * max(a, b)); = 2 ^ r; = p;

% Transformasi via FFT dengan zero padding Ff = fft2(Fs, p, q); Hf = fft2(H, p, q); % Konvolusi G = Hf .* Ff; % Peroleh citra hasil pemfilteran F = real(ifft2(G)); F = uint8(F(1:a, 1:b));

Akhir Program

Contoh penggunaan fungsi filterdft: >> H = fspecial(’sobel’);  >> F = filterdft('C:\Image\kotatua.png', H);  >> imshow(F)  6.8 Filter Lolos-Rendah Filter lolos-bawah (low-pass filter) adalah filter yang mempunyai sifat dapat meloloskan yang berfrekuensi rendah dan menghilangkan yang berfrekuensi tinggi.

Efek filter

ini membuat perubahan level keabuan menjadi lebih lembut.

Filter ini berguna untuk menghaluskan derau atau untuk kepentingan interpolasi tepi objek dalam citra. Tiga jenis filter lolos-rendah dilihat pada Tabel 6.1.


191

Tabel 6.1 Filter lolos-rendah Filter Ideal

Butterworth

Gaussian

Citra terkait

192


Jenis filter lolos-rendah pada kawasan frekuensi yang paling sederhana adalah yang dinamakan ILPF (Ideal Low Pass Filter). Filter ini memiliki fungsi transfer seperti berikut:

𝐻 (𝑣, 𝑢) = {

1 𝑗𝑖𝑘𝑎 𝐷(𝑣, 𝑢) ≤ 𝐷0 0 𝑗𝑖𝑘𝑎 𝐷(𝑣, 𝑢) > 𝐷0

(6.13)

Dalam hal ini, D0 adalah bilangan non-negatif yang biasa disebut radius filter, yang menentukan ambang frekuensi, dan D(v,u) adalah jarak antara (v,u) terhadap pusat filter, yang dinyatakan dengan 𝐷 (𝑣, 𝑢) = √𝑣 2 + 𝑢2

(6.14)

Skrip berikut menunjukkan suatu fungsi yang ditujukan untuk memfilter citra menggunakan ILPF. Program : ilpf.m

function F = ilpf(berkas, d0) % ILPF Digunakan untuk melaksanakan pemfilteran % pada kawasan frekuensi menggunakan ILPF. % Masukan: % berkas - nama citra % d0 - menentukan frekuensi ambang % Keluaran: % F - citra yang telah difilter Fs = double(imread(berkas)); [a, b] = size(Fs); %Peroleh ukuran citra % r p q


% Menentukan jangkauan frekuensi u dan v u = 0:(p - 1); v = 0:(q - 1); % Hitung indeks untuk meshgrid idx = find(u > q/2); u(idx) = u(idx) - q; idy = find(v > p/2); v(idy) = v(idy) - p;


193

% Peroleh array meshgrid [V, U] = meshgrid(v, u); % Hitung jarak D(v,u) D = sqrt(V.^2 + U.^2); % Hitung frekuensi ambang sebesar d0 kalai lebar citra ambang = d0 * q; % Peroleh fungsi transfer Hf = double(D <= ambang); % Transformasi via FFT dengan zero padding Ff = fft2(Fs, p, q); % Pemfilteran G = Hf .* Ff; % Transformasi balik F = real(ifft2(G)); F = uint8(F(1:a, 1:b));

Akhir Program

Contoh pemakaian fungsi ilpf: >> F = ilpf('C:\Image\kotatua.png', 0.08); imshow(F) 

Gambar 6.9 merupakan contoh penerapan IDLF untuk berbagai nilai d0 yang diterapkan pada kotatua.png.

194


Gambar 6.9 Hasil penerapan ILPF BLPF (Butterworth low pass filter) merupakan jenis filter lolos-rendah yang digunakan untuk memperbaiki efek bergelombang yang dikenal dengan sebutan ringing, yang diakibatkan oleh ILPF. Berbeda dengan ILPF, BLPF tidak memiliki titik diskontinu yang tajam. Fungsi transfernya berupa 𝐻 (𝑣, 𝑢) =

1 1+[𝐷 (𝑣,𝑢)/𝐷0 ]2𝑛

(6.15)

Dalam hal ini, n dinamakan orde filter. Skrip berikut menunjukkan suatu fungsi yang ditujukan untuk memfilter citra menggunakan BLPF. Program : blpf.m


function F = blpf(berkas, d0, n) % BLPF Digunakan untuk melaksanakan pemfilteran % pada kawasan frekuensi menggunakan BLPF. % Masukan: % berkas - nama citra % d0 - menentukan frekuensi ambang % n - menentukan faktor n % Keluaran: % F - citra yang telah difilter Fs = double(imread(berkas)); [a, b] = size(Fs); %Peroleh ukuran citra % r p q


% Menentukan jangkauan frekuensi u dan v u = 0:(p - 1); v = 0:(q - 1); % Hitung indeks untuk meshgrid idx = find(u > q/2); u(idx) = u(idx) - q; idy = find(v > p/2); v(idy) = v(idy) - p; % Peroleh array meshgrid [V, U] = meshgrid(v, u); % Hitung jarak D(v,u) D = sqrt(V.^2 + U.^2); % Menentukan n kalau n tidak disebutkan if nargin == 2 n = 1; end ambang = d0 * p; % Hitung frekuensi ambang Hf = 1 ./ (1 + D ./ ambang^(2 * n)); % Transformasi via FFT dengan zero padding Ff = fft2(Fs, p, q); % Pemfilteran G = Hf .* Ff; % Transformasi balik F = real(ifft2(G)); F = uint8(F(1:a, 1:b));

Akhir Program

195

196


Contoh pemakaian fungsi blpf: >> F = blpf('C:\Image\kotatua.png', 0.02, 0.3 );  >> imshow(F) 

Gambar 6.10 dan 6.11 merupakan contoh penerapan IDLF untuk berbagai nilai d0 dan n yang berbeda yang diterapkan pada kotatua.png.

Gambar 6.10 Hasil penerapan BLPF dengan N =1


197

Gambar 6.11 Hasil penerapan BLPF untuk berbagai N GLPF (Gaussian low pass filter) merupakan filter lolos-rendah dengan fungsi transfer seperti berikut: 𝐻 (𝑣, 𝑢) = 𝑒

𝐷2 (𝑣,𝑢) 2𝜎2

(6.16)

dengan 𝜎 merupakan deviasi standar. Sebagai contoh, dengan menggunakan 𝜎 sama dengan Do maka

𝐻 (𝑣, 𝑢) = 𝑒

𝐷2 (𝑣,𝑢) 2 2𝐷0

(6.17)

Saat D(v,u) = D0 , filter turun menjadi 0.607 terhadap nilai maksimum 1.

198


Contoh

berikut

merupakan

fungsi

yang

digunakan

pemfilteran dengan GLPF. Program : glpf.m

function F = glpf(berkas, d0) % GLPF Digunakan untuk melaksanakan pemfilteran % pada kawasan frekuensi menggunakan GLPF. % Masukan: % berkas - nama citra % d0 - menentukan frekuensi ambang % Keluaran: % F - citra yang telah difilter Fs = double(imread(berkas)); [a, b] = size(Fs); %Peroleh ukuran citra % r p q


% Menentukan jangkauan frekuensi u dan v u = 0:(p - 1); v = 0:(q - 1); % Hitung indeks untuk meshgrid idx = find(u > q/2); u(idx) = u(idx) - q; idy = find(v > p/2); v(idy) = v(idy) - p; % Peroleh array meshgrid [V, U] = meshgrid(v, u); % Hitung jarak D(v,u) D = sqrt(V.^2 + U.^2); % Menentukan n kalau n tidak disebutkan if nargin == 2 n = 1; end ambang = d0 * p; % Hitung frekuensi ambang Hf = exp(-(D.^2) ./ (2 * ambang ^ 2)); % Transformasi via FFT dengan zero padding Ff = fft2(Fs, p, q); % Pemfilteran G = Hf .* Ff; % Transformasi balik

untuk

melakukan


199

F = real(ifft2(G)); F = uint8(F(1:a, 1:b));

Akhir Program

Contoh pemakaian fungsi di glpf: >> F = glpf('C:\Image\kotatua.png', 0.05);  >> imshow(F) 

Gambar 6.12 merupakan contoh penerapan GLPF untuk berbagai nilai d0 yang diterapkan pada kotatua.png.

Gambar 6.12 Hasil penerapan filter GLPF

200


6.9 Filter Lolos-Tinggi Filter lolos-tinggi adalah filter yang ditujukan untuk menekan frekuensi rendah hingga frekuensi tertentu dan meloloskan frekuensi lainnya. Filter ini memiliki hubungan dengan filter lolos-rendah seperti berikut: 𝐻𝑙𝑡 (𝑣, 𝑢) = 1 − 𝐻𝑙𝑟 (𝑣, 𝑢)

(6.18)

Dengan Hlt (v,u) adalah fungsi transfer filter lolos-tinggi dan Hlf(v,u) adalah fungsi transfer filter lolos-rendah. Tiga jenis filter lolos-tinggi dilihat di Tabel 6.3. Ketiga filter yang tercantum dalam tersebut yaitu IHPF (Ideal high pass filter), BHPF (Butterworth high pass filter), dan GHPF (Gaussian high pass filter).

Tabel 6.3 Filter lolos-tinggi Filter Ideal

Butterworth

Citra terkait


Filter

201

Citra terkait

Gaussian

Contoh berikut menunjukkan implementasi filter BHPF untuk memfilter citra. Program : bhpf.m

function F = bhpf(berkas, d0, n) % BHPF Digunakan untuk melaksanakan pemfilteran % pada kawasan frekuensi menggunakan BHPF. % Masukan: % berkas - nama citra % d0 - menentukan frekuensi ambang % n - menentukan faktor n % Keluaran: % F - citra yang telah difilter Fs = double(imread(berkas)); [a, b] = size(Fs); %Peroleh ukuran citra % r p q


% Menentukan jangkauan frekuensi u dan v u = 0:(p - 1); v = 0:(q - 1); % Hitung indeks untuk meshgrid

202


idx = find(u > q/2); u(idx) = u(idx) - q; idy = find(v > p/2); v(idy) = v(idy) - p; % Peroleh array meshgrid [V, U] = meshgrid(v, u); % Hitung jarak D(v,u) D = sqrt(V.^2 + U.^2); % Menentukan n kalau n tidak disebutkan if nargin == 2 n = 1; end ambang = d0 * p; % Hitung frekuensi ambang Hlr = 1 ./ (1 + (D ./ ambang) .^(2 * n)); % Lolos-rendah Hlt = 1 - Hlr; % Lolos-tinggi % Transformasi via FFT dengan zero padding Ff = fft2(Fs, p, q); % Pemfilteran G = Hlt .* Ff; % Transformasi balik F = real(ifft2(G)); F = uint8(F(1:a, 1:b));

Akhir Program

Contoh pemakaian fungsi blpf: >> F = bhpf('C:\Image\goldhill.png', 0.005, 1 );  >> imshow(F) 

Gambar 6.13 menunjukkan gambar asli dan hasil pemrosesan dengan perintah di atas.


203

Gambar 6.13 Pemfilteran dengan BHPF Hasil pada Gambar 6.13(b) menunjukkan bahwa penerapan BHPF pada citra membuat latarbelakang menjadi hampir hilang, karena nilai intensitas reratanya hilang (menjadi nol). 6.10 Pemfilteran dengan Pendekatan High Frequency Emphasis Penerapan filter lolos-tinggi dengan cara yang telah dibahas menimbulkan efek berupa hilangnya latarbelakang. Hal ini disebabkan pemfilteran dengan cara tersebut menghilangkan komponen DC (F(0,0)). Nah,

untuk

mengatasi hal itu,

terdapat pendekatan yang dinamakan

pemfilteran high frequency emphasis (HFE). Dalam hal ini, penonjolan frekuensi tinggi diatur melalui rumus: 𝐻ℎ𝑓𝑒 (𝑣, 𝑢) = 𝑎 + 𝑏𝐻𝑙𝑡 (𝑣, 𝑢) Dalam hal ini, 

Hlt adalah fungsi transfer filter lolos-tinggi;



a adalah nilai ofset, sebagai penambah nilai rerata intensitas;



b adalah nilai pengali, untuk meningkatkan kontras.

(6.19)

204


Gonzales, dkk. (2004) menunjukkan bahwa penggunaan a sebesar 0,5 dan b sebesar 2 memberikan hasil yang memuaskan pada citra medis. Nilai a dan b itulah yang dicoba digunakan pada program berikut. Program : hfe.m

function F = hfe(berkas, d0, n) % HFE Digunakan untuk melaksanakan pemfilteran % pada kawasan frekuensi menggunakan BHPF % dan menerapkan HFE (High frequency emphasis). % Masukan: % berkas - nama citra % d0 - menentukan frekuensi ambang % n - menentukan faktor n % Keluaran: % F - citra yang telah difilter Fs = double(imread(berkas)); [a, b] = size(Fs); %Peroleh ukuran citra % r p q


% Menentukan jangkauan frekuensi u dan v u = 0:(p - 1); v = 0:(q - 1); % Hitung indeks untuk meshgrid idx = find(u > q/2); u(idx) = u(idx) - q; idy = find(v > p/2); v(idy) = v(idy) - p; % Peroleh array meshgrid [V, U] = meshgrid(v, u); % Hitung jarak D(v,u) D = sqrt(V.^2 + U.^2); % Menentukan n kalau n tidak disebutkan if nargin == 2 n = 1; end ambang = d0 * p; % Hitung frekuensi ambang Hlr = 1 ./ (1 + (D ./ ambang) .^(2 * n)); % Lolos-rendah Hlt = 1 - Hlr; % Lolos-tinggi % Proses HFE Hfe = 0.5 + 2 * Hlt;


205

% Transformasi via FFT dengan zero padding Ff = fft2(Fs, p, q); % Pemfilteran G = Hfe .* Ff; % Transformasi balik F = real(ifft2(G)); F = uint8(F(1:a, 1:b));

Akhir Program

Contoh penggunaan fungsi hfe: >> F = hfe('C:\Image\goldhill.png', 0.05, 1 );  >> imshow(F) 

Gambar 6.14 menunjukkan gambar asli goldhill.png dan hasil pemrosesan di atas.

Gambar 6.14 Hasil penerapan high frequency emphasis Terlihat bahwa terjadi penonjolan citra tanpa menghilangkan latarbelakang.

206


 Latihan 1. Berikan

beberapa

contoh

transformasi yang

mengalihkan ke

kawasan

frekuensi. 2. Terdapat data seperti berikut:

f(x) = (3, 4, 4, 5)

Hitunglah transformasi Fourier F(0) hingga F(3). 3. Jelaskan maksud istilah berikut. (a) Spektrum Fourier (b) Spektrum daya (power spectrum) 4. Apa yang Anda ketahui mengenai FFT? 5. Apa kegunaan fungsi fftshift? 6. Jelaskan yang dimaksud dengan istilah berikut. (a) Filter lolos-rendah (b) Filter lolos-tinggi 7. Apa yang dimaksud dengan wraparound error? Jelaskan cara mengatasinya. 8. Apa yang dimaksud dengan filter- filter berikut? (a) ILPF (b) BLPF (c) GLPF (d) IHPF (e) BHPF (f) GHPF (g) HFE 9. Buatlah fungsi bernama fillrb dengan 4 buah argumen yang secara berturut-turut menyatakan matriks citra, frekuensi ambang, dan orde filter. Fungsi tersebut memberikan nilai balik berupa fungsi transfer filter lolosrendah tipe Butterworth. Bantuan: modifikasilah berdasarkan blpf.m.


207

10. Gambarkan fungsi transfer fillrb pada Soal 9 dengan menggunakan f. Sebagai contoh, gunakan citra boneka2.png, dengan frekuensi sebesar 0,1, dan orde filter sebesar 1. Bantuan: Misalkan, H adalah fungsi transfer yang dihasilkan oleh fillrb, maka Anda bisa menyajikan filter transfer dengan menggunakan perintah berikut: mesh(fftshift(H)); colormap(jet); 11. Dengan menggunakan fungsi fillrb (Soal 9), bagaimana kalau dikehendaki untuk memperoleh fungsi transfer filter lolos-tinggi Butterworth?

208


Pengolahan Citra di Kawasan Frekuensi

Recommend Documents