Membuat Filter Digital dari Fungsi Transfer Dengan Transformasi Bilinear

admin

2 tahun ago

Pada tulisan ini akan diuraikan bagaimana membuat filter digital dalam bahasa C dari persamaan filter dalam transformasi Laplace. Contoh ini diadaptasi dari Example 8.4 di buku “Lizhe Tan, Digital Signal Processing, Fundamentals and Applications” halaman 324.

Misal diketahui suatu filter dalam domain s (Laplace) dengan persamaan sebagai berikut:

\(H(s)=\frac{10}{s+10} \)

Untuk mengubah ke bentuk digital, kita perlu tentukan perioda sampling. Misal dipilih perioda sampling T=0.01 detik

Persamaan transformasi s ke z sebagai berikut

\( s=\frac{2(z -1) }{Tz + 1} \)

H(s) diubah ke H(z) dengan cara substitusi s pada H(s). T=0.01 dimasukkan ke persamaan transformasi. Maka didapat persamaan berikut:

\( H(z)=\frac{10}{\frac{200(z-1)}{z+1}+10 } \)

\( H(z)=\frac{0.05}{\frac{z-1}{z+1}+0.05 } \)

\( H(z)=\frac{0.05(z+1)}{z-1+0.05(z+1) } \)

\( H(z)=\frac{0.05z+0.05}{1.05z – 0.95 } \)

\( H(z)=\frac{(0.05z+0.05)/1.05z}{(1.05z – 0.95)/1.05z } \)

Didapatkan H(z) dengan bentuk standar.

\( H(z)=\frac{0.0476+0.0476z^{-1}}{1 – 0.9048 z^{-1} } \)

Selanjutnya terapkan teknik di Bab 6 (Example 6.5) untuk mengubah persamaan dalam domain (z) ke domain waktu (n).

Akan didapat persamaan difference dalam n sebagai berikut:

y(n)=0.0476 x(n) + 0.0476 x(n-1) + 0.9048 y(n-1)

Simulasi

Tahap selanjutnya adalah simulasi persamaan difference tersebut dengan input fungsi step.

Tahap simulasi ini penting untuk mengecek apakah filter digital yang dihasilkan memiliki perilaku yang sama dengan filter analog di awal. Jika hasilnya berbeda, ada kemungkinan salah perhitungan atau salah koding.

Simulasi dapat dilakukan dengan bahasa pemrograman apa saja, misal C, C++, Matlab, Python, dan sebagainya.

Simulasi di bawah ini dibuat dengan bahasa Python. GUI yang dipakai adalah Jupyter Notebook

import matplotlib
import matplotlib.pyplot as plt
import numpy as np
from matplotlib.ticker import AutoMinorLocator
from matplotlib import gridspec
import matplotlib.ticker as ticker

%matplotlib inline

t=-1 # start time at -1
xn=0
xn1=0
yn1=0
period=0.01

times=np.linspace(-1, 10, 111).reshape(111, 1)
xs=[]
ys=[]
ts=[]
for counter in range(0,200):
    if t<0:
        xn=0        
    else:
        xn=1    
    t=t+period
    xs.append(xn)                    
    y=0.04676 * xn + 0.0476 * xn1 + 0.9048 * yn1
    ys.append(y)
    ts.append(t)
    yn1=y
    xn1=xn

fig = plt.figure(figsize=(8,4))
gs = gridspec.GridSpec(1,1)
ax1 = fig.add_subplot(gs[0])
ax1.plot(ts,xs,"g-",label="input x")
ax1.plot(ts,ys,"r-",label="output y")
ax1.plot(1/10,.632,"bo",label="output y")
ax1.set_xlabel("waktu",family="serif",  fontsize=12)
ax1.set_ylabel("value",family="serif",  fontsize=12)
ax1.legend(loc='best')
fig.tight_layout()
plt.grid(True)
fig.savefig("blt-step.jpg", format="jpg",dpi=65)

Berikut ini output sistem (merah) terhadap input step (warna hijau).

time constant = 1/10

Titik biru adalah t = time constant, dengan value = 63.2%. Posisi titik biru ini pada grafik merah, sehingga dapat disimpulkan respon sistem ini cocok dengan yang diharapkan.

Selanjutnya adalah memberi sinyal sinusoidal dengan frekuensi = frekuensi cut-off ke filter tersebut . Menurut teori, amplitudo pada steady state adalah 0.707 x amplitudo maksimum.

Berikut ini kode untuk memberikan sinyal sinusoida.

t=-1 # start time at -1
xn=0
xn1=0
yn1=0
period=0.01

times=np.linspace(-1, 10, 111).reshape(111, 1)
xs=[]
ys=[]
ts=[]
w=10
for counter in range(0,200):    
    t=t+period
    xn=math.sin(w*t)
    xs.append(xn)                    
    y=0.04676 * xn + 0.0476 * xn1 + 0.9048 * yn1
    ys.append(y)
    ts.append(t)
    yn1=y
    xn1=xn

Berikut ini kode untuk menampilkan grafik sinyal

fig = plt.figure(figsize=(8,4))
gs = gridspec.GridSpec(1,1)
ax1 = fig.add_subplot(gs[0])
ax1.plot(ts,xs,"g-",label="input x")
ax1.plot(ts,ys,"r-",label="output y")
ax1.set_xlabel("waktu",family="serif",  fontsize=12)
ax1.set_ylabel("value",family="serif",  fontsize=12)
ax1.legend(loc='best')
fig.tight_layout()
plt.grid(True)
fig.savefig("blt-sinusoidal.jpg", format="jpg",dpi=65)

Grafik tersebut menunjukkan amplitudo output (merah) ada di sekitar 0.7 dari amplitudo sinyal input (hijau). Terjadi pergeseran fasa sekitar 45 derajat.

Implementasi Bahasa C di Windows

Tahap selanjutnya adalah membuat implementasi filter dengan bahasa C.

Compiler yang digunakan adalah CodeBlocks

Berikut ini kode filter digital dalam bahasa C.

Output program ini adalah file win10-filter.csv

#include <stdio.h>
#include <stdlib.h>

// inisialisasi filter digital

void filter_init(float *xn1,float *yn1)
{
    *yn1=0;
    *xn1=0;
}

// hitung output filter digital filter digital
// delayed variabel diberikan sebagai parameter
float filter_run(float *xn,float *xn1,float *yn1)
{
    float y;
    y=0.04676 * *xn + 0.0476 * *xn1 + 0.9048 * *yn1;
    *yn1=y;
    *xn1=*xn;
    return y;
}

int main()
{
    FILE * fp;
    float y,xn,xn1,yn1;
    int i;
    float t;
    float period;
    t=-1;
    period=0.01;

    fp = fopen ("win10-filter.csv", "w+");
    filter_init(&xn1,&yn1);

    for (i=0; i<200; i++)
    {
        t=t+period;
        if(t<0)
        {
            xn=0;
        }
        else
        {
            xn=1;
        }
        y=filter_run(&xn,&xn1,&yn1);
        fprintf(fp,"%f,%f,%f\n",t,xn,y);
    }
    fclose(fp);
    return 0;
}

Berikut ini grafik dari file csv dengan Excel

Bentuk grafik sesuai dengan grafik di Python. Nilai output csv dapat dibandingkan dengan simulasi di Python. Hasilnya sesuai.

Implementasi Bahasa C di ATmega328

Implementasi filter di ATmega328 diperlukan jika kita membuat filter digital atau membuat simulator plant untuk Hardware In the Loop (HIL)

under construction

Implementasi Bahasa C di ESP32

Implementasi filter di ESP32diperlukan jika kita membuat filter digital atau membuat simulator plant untuk Hardware In the Loop (HIL)