Dalam suatu sistem perekaman suara, dilakukan perekaman suara secara digital dalam suatu area yang berisik dengan frekuensi cuplik 8 kHz. Dianggap bahwa rekaman suara yang dihasilkan mengandung informasi mulai dari hingga 1800 Hz, sehingga kita bisa merancang sebuah tapis lolos-rendah (low-pass) yang akan menahan derau antara 1800 Hz hingga batasan Nyquist, 4000 Hz (separo dari frekuensi cuplik 8000 Hz).
Berikut spesifikasi Tapis kita:
- Tipe Tapis = FIR lolos-rendah (lowpass FIR filter)
- Jangkauan Frekuensi Passband = 0 – 1800 Hz
- Riak (ripple) pada Passband = 0,02 dB
- Jangkauan Frekuensi Stopband = 2000 – 4000 Hz
- Pelemahan (attenuation) pada Stopband = 50 dB
- Tipe Jendela = Hamming
- Jumlah koefisien = 133
- Frekuensi cutoff = 1900 Hz
Fungsi jendela menggunakan persamaan-persamaan berikut:
Selanjutnya dibuat skrip fungsi firwd(N,Ftype,WnL,WnH,Wtype) yang digunakan untuk membuat tapis menggunakan jendela dengan listing sebagai berikut…
function B=firwd(N,Ftype,WnL,WnH,Wtype)
% B = firwd(N,Ftype,WnL,WnH,Wtype)
% Rancang Tapis FIR menggunakan metode Fungsi Jendela
% Parameter2 masukan:
% N: jumlah koefisien (tap) tapis FIR
% Catatan: Harus bilangan ganjil.
% Ftype: Tipe tapis
% 1. Tapis Low-pass ;
% 2. Tapis High-pass;
% 3. Tapis Bandpass;
% 4. Tapis Bandstop.
% WnL: Frekuensi cutoff bawah (radian), WnL=0 untuk Tapis Highpass
% WnH: Frekuensi cutoff atas (radian), WnH=0 untuk Tapis Lowpass.
% Wtype: tipe fungsi jendela
% 1. Jendela Rectangular;
% 2. Jendela Triangular;
% 3. Jendela Hanning;
% 4. Jendela Hamming;
% 5. Jendela Blackman;
% Keluaran:
% B: Koefisien2 Tapis FIR.
%
M=(N-1)/2;
hH=sin(WnH*[-M:1:-1])./([-M:1:-1]*pi);
hH(M+1)=WnH/pi;
hH(M+2:1:N)=hH(M:-1:1);
hL=sin(WnL*[-M:1:-1])./([-M:1:-1]*pi);
hL(M+1)=WnL/pi;
hL(M+2:1:N)=hL(M:-1:1);
if Ftype == 1
h(1:N)=hL(1:N);
end
if Ftype == 2
h(1:N)=-hH(1:N);
h(M+1)=1+h(M+1);
end
if Ftype ==3
h(1:N)=hH(1:N)-hL(1:N);
end
if Ftype == 4
h(1:N)=hL(1:N)-hH(1:N);
h(M+1)=1+h(M+1);
end
% Fungsi-fungsi jendela;
if Wtype ==1
w(1:N)=ones(1,N);
end
if Wtype ==2
w=1-abs([-M:1:M])/M;
end
if Wtype ==3
w= 0.5+0.5*cos([-M:1:M]*pi/M);
end
if Wtype ==4
w=0.54+0.46*cos([-M:1:M]*pi/M);
end
if Wtype ==5
w=0.42+0.5*cos([-M:1:M]*pi/M)+0.08*cos(2*[-M:1:M]*pi/M);
end
B=h .* w
Fungsi firwd() ini akan kita gunakan dalam ekperimen Noise Removal kita ini. Baik, mari kita lakukan proses tapis derau langkah demi langkah sebagai berikut, kita awali dengan menutup jendela-jendela, sekalian membersihkan memori Matlab…
close all; clear all
Kita siapkan frekuensi cuplik (8000 Hz) dan periodenya (1/Fs), sekaligus membaca data yang digunakan untuk eksperimen ini (menggunakan instruksi load pada Matlab, data beserta listing program lengkap bisa diunduh disini):
fs=8000;T=1/fs;
load we.dat
Kemudian kita siapkan variabel waktunya (t). Mulai dari 0 (detik) hingga [Panjang(we)-1] dikali T (detik)…
t=[0:length(we)-1]*T;
Baik, sekarang kita siapkan deraunya. Derau berkaitan dengan data-data, tingkatnya tinggi (high level) dan bersifat broadband (nanti akan Anda lihat spektrumnya), diawali dengan 1/4 dari rerata kuadrat data…
th=mean(we.*we)/4;
v=sqrt(th)*randn([1,length(we)]);
Gabungkan derau (v) dengan sinyal aslinya (we) dan gambarkan hasilnya (dalam ranah waktu)….
x=we+v;
subplot(2,1,1);plot(t,x,’k’);
xlabel(‘Waktu (detik)’);ylabel(‘Amplitudo’);grid;
Kemudian gambarkan juga spektrumnya…
N=length(x);
f=[0:N/2]*fs/N;
Axk=2*abs(fft(x))/N;Axk(1)=Axk(1)/2;
subplot(2,1,2); plot(f,Axk(1:N/2+1),’k’);
xlabel(‘Frekuensi (Hz)’); ylabel(‘Besaran |X(f)| ‘);grid;
hasilnya…
Bagaimana? Sekarang Anda bisa melihat dengan jelas sinyal + derau dan spektrumnya. Perhatikan bentuk-bentuk ‘rumput’ pada spektrumnya (ini karena adanya derau). Selanjutnya dilakukan proses penapisan untuk mengambil (atau lebih tepatnya menahan) komponen derau. Kita siapkan frekuensi cutoff dalam radian yang dinormalisasi terhadap frekuensi cuplik (Fc/Fs) dikali dengan 2*pi…
wc=2*pi*1900/8000;
Kemudian kita siapkan tapisnya…
B=firwd(133,1,wc,0,4);
Dan akhirnya kita lakukan proses penapisan (filtering)…
y=filter(B,1,x);
Selanjutnya kita gambarkan hasilnya, sekalian spektrumnya…
Ayk=2*abs(fft(y))/N;Ayk(1)=Ayk(1)/2;
subplot(2,1,1); plot(t,y,’k’);
xlabel(‘Waktu (detik)’);ylabel(‘Amplitudo’);grid;
subplot(2,1,2);plot(f,Ayk(1:N/2+1),’k’);
xlabel(‘Frekuensi (Hz)’); ylabel(‘Besaran |Y(f)| ‘);grid;
Hasilnya…
Okey! Bandingkan kedua gambar tersebut, adakah perbedaan? Dalam ranah waktu? Dalam ranah frekuensi? Mo Jawab? Komentar? Pertanyaan? Silahkan gunakan form komentar di bawah ini…
Terima kasih, semoga bermanfaat.. Oya mo listingnya lengkap (termasuk data we.dat)? Unduh aja disini (RAR)!
Catatan:
- Beberapa instruksi yang tidak dijelaskan dalam artikel mohon merujuk HELP-nya Matlab.
sumber:
- Tan, Li, 2008, Digital Signal Processing: Fundamentals and Applications, Academic Press, Elsevier, USA