TKE 2403
SISTEM PENGOLAHAN ISYARAT
Kuliah 8 – Analisis Wavelet : Wavelet Kontinyu
Indah Susilawati, S.T., M.Eng.
Program Studi Teknik Elektro Fakultas Teknik dan Ilmu Komputer Universitas Mercu Buana Yogyakarta 2009
KULIAH 8 SISTEM PENGOLAHAN ISYARAT
Analisis Wavelet : Wavelet Kontinyu Pendahuluan Sebagaimana telah diuraikan pada bab sebelumnya, sebuah alihragam dapat dipandang sebagai suatu pemetaan kembali suatu sinyal tertentu sehingga dapat memperlihatkan lebih banyak informasi daripada sinyal aslinya. Alihragam Fourier memberikan informasi mengenai konten frekuensi yang sebelumnya tidak tampak jika sinyal disajikan dalam domain waktu (sinyal asli). Namun alihragam Fourier tidak mampu mendiskripsikan karakteristik waktu dan frekuensi secara bersamaan. Alihragam wavelet dapat digunakan untuk mendiskripsikan
suatu sinyal yang ubah waktu; sinyal tidak dibagi menjadi
bagian-bagian menurut waktu namun dibagi menjadi segmen-segmen menurut suatu skala. Alihragam Fourier dinyatakan sebagai
X (ω m ) =
∞
− jωm t x ( t ) e dt ∫
(1)
−∞
Sinyal sinusoidal biasanya digunakan dalam aliragam Fourier karena sifat frekuensinya yang unik yaitu energinya terkandung pada satu frekuensi tertentu. Sebenarnya sebarang fungsi (selain sinusoidal) dapat digunakan. Jika fungsi tersebut mempunyai durasi yang terhingga saja, maka akan dapat digeser-geser sepanjang sinyal x(t) (seperti dilakukan pada proses konvolusi) dan dihasilkan alihragam Fourier yang sifatnya short-term atau STFT (short-term Fourier transform). ∞
STFT (t , f ) = ∫ x(τ ) ( w(t − τ ) e − j 2π f τ ) dτ −∞
(2)
Dengan f adalah frekuensi yang juga mengindikasikan keanggotaan, w(t - τ) adalah fungsi jendela dengan t berlaku sebagai penggeser jendela sepanjang x. Secara umum persamaan (2) dapat ditulis kembali sebagai persamaan berikut. ∞
X (t , m) = ∫ x(τ ) f (t − τ ) m dτ
(3)
−∞
Dengan f(t)m adalah fungsi keanggotaan dan m menspesifikasikan nomor keanggotaan.
Alihragam Wavelet Kontinyu Alihragam wavelet mempergunakan konsep yang dinyatakan pada persamaan (3). Dalam analisis wavelet, ada sejumlah fungsi yang dapat digunakan, namun keanggotaannya selalu terdiri atas versi-versi yang diregangkan atau dimampatkan dari fungsi basisnya, sebagaimana halnya pergeseran. Konsep ini mengarah pada definisi persamaan untuk alihragam wavelet kontinyu (CWT: Continuous Wavelet Transform), ∞
W ( a, b) =
∫ x(τ )
−∞
1 ⎛t −b⎞ ψ *⎜ ⎟ dt a |a| ⎝ ⎠
(4)
Dengan b berlaku sebagai penggeser fungsi sepanjang x(t) dan a berlaku sebagai penskala waktu pada fungsi ψ. Jika a lebih besar daripada 1, fungsi wavelet ψ menjadi terregang sepanjang sumbu waktu dan jika a kurang daripada 1 (dan masih bernilai positif), maka ψ akan termampatkan. Jika nilai a negatif maka akan membuat fungsi wavelet ψ menjadi terbalik pada arah sumbu waktu. Tanda * menunjukkan operasi konjugat kompleks, dan normalisasi
1 |a|
memastikan bahwa energinya sama untuk semau nilai a (dan juga untuk semua nilai b). Jika b = 0 dan a = 1, maka wavelet berada dalam bentuk naturalnya dan disebut dengan istilah the mother wavelet yaitu ψ1,0(t) = ψ(t). Gambar 1 memperlihatkan hal ini (dan beberapa anggotanya yang dihasilkan oleh
beberapa pergeseran yaitu dilasi dan kontraksi). Pada gambar tersebut digunakan wavelet Morlet yang dinyatakan dengan
⎛
2
⎞
ψ (t ) = e −t cos ⎜⎜ π t ⎟⎟ ln 2 ⎝ ⎠ 2
(5)
Gambar 1. Mother wavelet (a=1) dengan dua dilasi (a = 2 dan 4) dan satu kontraksi (a = 0,5)
Koefisien wavelet, W(a, b), menggambarkan korelasi antara sinyal dengan wavelet pada berbagai pergeseran dan skala: kemiripan antara sinyal dan wavelet pada suatu kombinasi skala dan posisi yaitu a,b. Dengan kata lain, koefisien wavelet menentukan amplitude deret wavelet pada jangkauan skala dan pergeseran yang harus dijumlahkan untuk merekontruksi kembali sinyal. Jika fungsi wavelet ψ(t) dapat dipilih dengan tepat, maka dimungkinkan untuk merekontruksi sinyal aslinya dari koefisien wavelet (seperti yang dilakukan pada alihragam Fourier). Oleh karena CWT mendekomposisi sinyal menjadi koefisien-
koefisien dari dua variabel yaitu a dan b, maka dibutuhkan penjumlahan (atau integrasi) ganda untuk memperoleh kembali sinyal asli dari koefisien-koefisien tersebut. Hal ini dinyatakan pada persamaan berikut. ∞
1 x(t ) = C a =∫−∞
∞
∫ W ( a, b) ψ
a ,b
(t ) da db
(6)
b = −∞
Dengan ∞
| ψ (ω ) | 2 C= ∫ dω |ω | −∞ dan 0 < C < ∞
(disebut kondisi yang diijinkan atau admissibility condition).
Biasanya jika dikehendaki untuk memperoleh kembali sinyal asli maka akan digunakan alihragam wavelet diskret (karena lebih mudah dan perhitungan yang harus dilakukan juga lebih sedikit).
Karakteristik Waktu – Frekuensi Wavelet Wavelet seperti yang ditunjukkan pada Gambar 1 tidak berada pada satu waktu atau frekuensi yang spesifik. Kenyataannya wavelet menawarkan kompromi antara letak waktu dan frekuensi: wavelet terlokalisir dalam waktu dan frekuensi (namun tidak untuk salah satunya). Ukuran jangkauan waktu untuk suatu wavelet tertentu, yaitu Δtψ, dapat ditentukan oleh akar kuadrat momen kedua (second moment) wavelet tersebut terhadap pusat waktunya (yaitu momen pertama),
(7) Dengan t0 adalah waktu pusat atau momen pertama dari wavelet dan ditentukan dengan persamaan berikut.
(8) Dengan cara yang sama maka jangkauan frekuensi suatu wavelet, yaitu Δωψ, dapat ditentukan dengan
(9) Dengan ψ(ω) adalah representasi domain frekuensi (yaitu alihragam Fourier) dari ψ(t/a) dan ω0 adalah frekuensi pusat dari ψ(ω). Frekuensi pusat ditentukan dengan persamaan yang serupa dengan persamaan (8), yaitu
(10) Jangkauan waktu dan frekuensi dari keanggotaan dapat diperoleh dari mother wavelet menggunakan persamaan (7) dan (9). Dilasi oleh variabel a mengubah jangkauan waktu dengan cara mengalikan Δtψ dengan a. Dengan demikian, jangkauan waktu dari ψa,0 didefinisikan sebagai Δtψ (a) = |a| Δtψ. Hubungan terbalik antara waktu dan frekuensi diperlihatkan pada Gambar 2 yang diperoleh dengan menerapkan persamaan (7 – 10) pada wavelet Mexican Hat. Wavelet Mexican Hat dinyatakan dengan persamaan berikut.
(11)
Gambar 2. Batas-batas waktu-frekuensi dari wavelet Mexican Hat untuk berbagai nilai a Jangkauan frekuensi, atau bandwidth, akan sama dengan jangkauan mother wavelet dibagi a; yaitu Δωψ (a) = Δωψ / |a|. Jika jangkauan frekuensi dikalikan dengan a, maka akan dihasilkan sebuah konstanta yaitu hasil kali konstanta yang dihasilkan oleh persamaan (7) dan (9), (12) Persamaan (12) menunjukkan bahwa hasil kali tersebut invarian terhadap dilasi dan bahwa jangkauan terhubung terbalik: kenaikan pada jangkauan frekuensi, Δωψ (a) akan menurunkan jangkauan waktu Δtψ (a). Jangkauanjangkauan ini berhubungan dengan resolusi waktu dan frekuensi dari CWT. Dengan menurunkan jangkauan waktu wavelet (yaitu dengan menurunkan a) akan menghasilkan karakteristik waktu yang lebih akurat namun dengan mengorbankan resolusi frekuensi. Demikian juga jika sebaliknya.
Oleh karena resolusi waktu dan frekuensi mempunyai hubungan terbalik, maka CWT akan menawarkan resolusi yang lebih baik jika a besar dan panjang wavelet (dan jendela waktu efektifnya) panjang. Sebaliknya, jika a kecil, waveletnya pendek dan resolusi waktu-nya maksimum namun wavelet hanya merespon komponen frekuensi tinggi saja. Karena a bervariasi, maka dapat dilakukan tadeoff antara resolusi waktu dan frekuensi, dan ini merupakan kunci CWT.
Contoh 1 Tulis program dalam Matlab untuk membentuk CWT dari sebuah sinyal yang terdiri atas dua sinyal sinus yang berturutan dengan frekuensi 10 dan 40 Hz (seperti terlihat pada Gambar 3). Plot koefisien wavelet sebagai fungsi a dan b. Gunakan wavelet Morlet.
Sinyal Sinus 10 dan 40 Hz 1 0.8 0.6 0.4
amplitude
0.2 0 -0.2 -0.4 -0.6 -0.8 -1
0
0.5
1
1.5
2
2.5
waktu, detik
Gambar 3. Dua sinyal sinus berbeda frekuensi (10 dan 40 Hz) dengan perubahan frekuensi secara mendadak (step)
% Contoh 1 dan Gambar 3 dan Gambar 4 % Membangkitkan dua sinusoaida yang berubah frekuensi secara step % Menerapkan alihragam wavelet kontinyu dan menggambar hasilnya clear all; close all; % Setting konstanta fs = 500; N = 1024; N1 = 512; f1 = 10; Hz f2 = 40; Hz resol_level = 128; decr_a = .5; a_init = 4; wo = pi*sqrt(2/log2(2));
% % % %
frekuensi sampling panjang sinyal jumlah titik wavelet frekuensi sinusoida pertama dalam
% frekuensi sinusoida kedua dalam % % % %
jumlah nilai a penurunan a nilai awal a faktor skala frekuensi wavelet
% Membangkitkan dua gelombang sinus pada Gambar 3 tn = (1:N/4)/fs; % vektor waktu utk membangkitkan sinus b = (1:N)/fs; % vektor waktu untuk plotting x = [zeros(N/4,1);sin(2*pi*f1*tn)';sin(2*pi*f2*tn)';zeros(N/4,1)]; plot (b,x) xlabel('waktu, detik') ylabel('amplitude') title('Sinyal Sinus 10 dan 40 Hz') ti = ((1:N1/2)/fs)*10; % vektor waktu untuk membangun wavelet % Menghitung alihragam wavelet kontinyu (wavelet morlet) for i = 1:resol_level a(i) = a_init/(1+i*decr_a); % set skala t = abs(ti/a(i)); % vektor skala untuk wavelet mor = (exp(-t.^2).*cos(wo*t))/sqrt(a(i)); wavelet = [fliplr(mor) mor] % membuat simetri terhadap nol ip = conv(x, wavelet); % konvolusi wavelet dan sinyal ex = fix((length(ip)-N)/2); % menghitung titik tambahan/2 CW_Trans(:,i) = ip(ex+1:N+ex,1); end % Plot dalam 3 dimensi figure d = fliplr(CW_Trans); mesh(a,b,CW_Trans) xlabel('skala a') ylabel('b (detik)') zlabel('CWT') title('CWT Sinyal Sinus 10 dan 40 Hz') rotate3d on % view(130,50)
Ket: Nama file SPI_wavelet_1.m
Gambar 4a. Koefisien wavelet dari CWT sinyal yang terdiri atas dua sinyal sinus sekuensial dengan frekuensi 10 dan 40 Hz menggunakan wavelet Morlet
Gambar 4b. Gambar 4a yang pada sudut azimuth (az) 130 dan elevasi (el) 50
Contoh 2 Tentukan batas-batas waktu – frekuensi dari wavelet Mexican Hat dengan menerapkan persamaan (7 – 10). Untuk setiap nilai a, wavelet yang terskala dibentuk menggunakan pendekatan yang sama dengan yang digunakan pada contoh 1. Magnitude kuadrat dari tanggapan frekuensi dihitung menggunakan FFT. Waktu pusat, t0, dan frekuensi pusat, w0, dibentuk dengan menggunakan persamaan (8) dan (10) secara langsung. Perlu dicatat bahwa karena wavelet dibentuk simetris terhadap t = 0 maka waktu pusat t0 akan selalu nol dan waktu offset t1 ditambahkan untuk kepentingan plotting. Batas-batas waktu dan frekuensi dihitung menggunakan persamaan (7) dan (9), hasilnya digambar sebagai sebuah persegipanjang yang bersesuaian dengan waktu pusat dan frekuensi pusat masing-masing.
% % % %
Contoh 2 dan Gambar 5 Menggambar batas-batas wavelet untuk berbagai nilai 'a' Menentukan jangkauan waktu dan frekuensi dari wavelet Mexican Hat Menggunakan persamaan (7 - 10)
clear all; close all; N = 1000; fs = 1000; (asumsi) wo1 = pi*sqrt(2/log2(2)); wavelet a = [.5 1.0 2.0 3.0]; x1 = ((1:N/2)/fs)*10; t = (1:N)/fs; omega = (1:N/2)*fs/N; for i = 1:length(a) t1 = x1./a(i) mex = exp(-t1.^2).*(1-2*t1.^2); w = [fliplr(mex) mex]; wsq = abs(w).^2 W = fft(w); frekuensi dan Wsq = abs(W(1:N/2)).^2; digunakan
% Panjang data % Frekuensi pencuplikan % Konstanta utk skala waktu % % % %
Nilai-nilai a Wavelet selama +/- 10 detik Skala waktu Skala frekuensi
% % % % %
Vektor waktu utk wavelet Wavelet Mexican Hat Simetri terhadap nol Kuadrat dari wavelet Mencari representasi
% kuadratnya. Jangkauan yg
% hanya fs/2 t0 = sum(t.*wsq)/sum(wsq); % Menghitung waktu pusat d_t = sqrt(sum((t-t0).^2.*wsq)/sum(wsq)); % Menghitung sebaran waktu w0 = sum(omega.*Wsq)/sum(Wsq); % Menghitung frekuensi pusat d_w0 = sqrt(sum((omega-w0).^2.*Wsq)/sum(Wsq));
t1 = t0*a(i);
% Mengatur posisi waktu
hold on; % Menggambar batas-batas waktu - frekuensi plot([t1-d_t t1-d_t], [w0-d_w0 w0+d_w0],'m'); plot([t1+d_t t1+d_t], [w0-d_w0 w0+d_w0],'m'); plot([t1-d_t t1+d_t], [w0-d_w0 w0-d_w0],'m'); plot([t1-d_t t1+d_t], [w0+d_w0 w0+d_w0],'m'); title('Batas-batas waktu-frekuensi wavelet Mexican Hat') xlabel('waktu (detik)') ylabel('frekuensi (rad/det)') end
Ket: Nama file SPI_wavelet_2.m
Batas-batas waktu-frekuensi wavelet Mexican Hat 11 10 9
frekuensi (rad/det)
8 7 6 5 4 3 2 1 0.2
0.4
0.6
0.8
1 1.2 waktu (detik)
1.4
1.6
1.8
Gambar 5. Batas-batas waktu-frekuensi dari wavelet Mexican Hat untuk berbagai nilai a (bandingkan dengan yang diperlihatkan pada Gambar 2).