淺談動平均濾波器

動平均濾波器 (moving average filter) 可說是世界上使用最廣泛的一種濾波器了; 在股市裡面, 每天都可以看到動平均濾波器的應用, 股價趨勢圖裡面的週線、月線與年線, 都是長度不等的動平均輸出圖形。我們在實作訊號處理的時候, 也很常應用這個簡單又有效率的動平均濾波器。 由於它應用廣泛, 許多程式設計師也許沒有修過數位訊號處理的課程, 不熟悉數位訊號處理分析方式, 但是又需要寫這樣的程式。因此, 在這一篇文章中, 筆者將會儘量以實例來介紹動平均濾波器的概念, 幾種不同的實作方式以及優缺點。

詳全文請點: 全文 PDF

定點數表示法

筆者今天要來介紹定點數表示法, 也就是 fixed-point representation. 讀者可能會很納悶, 現在有浮點數的機器非常普遍, 而且速度都很快, 比如說 Intel Core i7 處理器,開啟 AVX SIMD 的情況下, 可以一次處理四個 64 位元浮點數的乘法。但是在很多地方, 比如說訊號處理、控制系統、或者很多手機內部的處理器, 都欠缺浮點運算能力,這時候就必須使用定點數。...全文 PDF

用 Linear Filter(IIR) 來實作 Fibonacci 序列

PyCon.TW 2012 在日前圓滿落幕, 會中的 keynote speaker Travis Oliphant 講
了一個使用 linear filter 實作 Fibonacci 的例子, 在演講中只有一張投影片,
在這裡我嘗試做分解動作解釋一下。投影片內容在 http://www.slideshare.net/pycontw/largescale-arrayoriented-computing-with-python

我們知道 Fibonacci 序列是 x(n) = x(n-1) + x(n-2), 而 x(0) = 0, x(1)
= 1. Fibonacci 序列的產生也常常是教科書介紹 recursive 的好範例, 不過就
如 Travis Oliphant 所示, 用 recursive 實作的複雜度是 exponential 成長的,
而我們也可以用 SciPy 裡面的 linear filter 來實作這件事, 方法如下：

from scipy.signal import lfilter
from numpy import zeros
b = array([1.0])
a = array([1., -1, -1])
zi = array([0, 1.0])
def fib3a(N):
   y, zf = lfilter(b, a, zeros(N, dtype=float), zi=zi)
   return y

我來慢動作分解一下上面的實作方法。首先 lfilter 在 b = array([1.0]) 以及
a=array([1.,-1,-1]) 的情況下會產生b(0)x(n) = a(0)y(n) + a(1)y(n-1) +
a(2)y(n-2) 的式子, 代入 a,b, 我們可以得到 x(n) = y(n) - y(n-1) -
y(n-2),而 lfilter 的第三個參數就是 x 的序列, 是 zeros(N, dtype=float),
也就是說輸入全為零,上述的式子就會變成 0 = y(n) - y(n-1) - y(n-2), 把 y(n) 移到等式的左方, 就可以得到 y(n) = y(n-1) + y(n-2), 這就是 Fibonacci 數列的表示法啦! 最後的 y(n) 就是 Fibonacci 的第 n 個數列。另外,
zi=array([0,1.0]) 就是當 lfilter 開始執行的時候, delay element 裡面的元
素, 也就是 y(0) 與 y(1) 的值。

以濾波器的角度而言, 這個濾波器在輸入全為零的狀況下, 自己產生
0,1,1,2,3,5,… 這種發散的數列, 等於自己在震盪, 是我們在實際上不會去使
用的濾波器, 拿來實作 Fibonacci 數列, 實在是有趣。

--

LaTeX PSTricks 訊號處理範例: z-Transform

上一篇的 moving average filter y(n) = (15/16)*y(n-1) + (1/16)*x(n) 的
z-transform 有一個 pole 在 z=(15/16, 0) 的點上, 要畫這個 z-Transform 在
z-plane 上的表示圖, 可以用下列的 PSTricks 來表示。

\usepackage{pstricks}
\usepackage{pst-sigsys}

圖形如下表示:

\begin{center}
  \begin{pspicture}[showgrid](-2,-2)(2,2)
    \pscircle[linecolor=gray](0,0){1} % unit circle
    \pspole(0.9375,0){z1}
    \nput{0}{z1}{$(\frac{15}{16},0)$}
  \end{pspicture}
\end{center}

--

LaTeX PSTricks 訊號處理範例: Moving Average Filter

今天要介紹的是如何用 PSTricks 來繪製訊號處理的 functional diagram, 假設
我們要描述的訊號處理方程式是: y(n) = (15/16)*y(n-1) + (1/16)*x(n). 這是
一個動平均濾波器, y(n) 是目前的輸出值, x(n) 是目前的輸入值, y(n-1) 則是上一次的輸出值,
那他的 functional diagram 會是如下圖所示:

上述圖形使用了 PSTricks 以及 PSTricks 的 pst-sigsys 套件，所以要在 preamble 的部份如下宣告:

\usepackage{pstricks}
\usepackage{pst-sigsys}

描述這個濾波器的 LaTeX source code 如下:

\begin{center}
  % To describe y(n) = (15/16)*y(n-1) + (1/16)*x(n)
  \begin{pspicture}(-2,-1)(10,2)
    \pssignal(0,1){x}{$x(n)$}
    \pscircleop[oplength=0.25,operation=times](2,1){op1}
    \pssignal(2,0){coefx}{$\frac{1}{16}$}
    \pscircleop[oplength=0.25](4,1){op2}
    \pssignal(8,1){y}{$y(n)$}
    \dotnode(6,1){ydot}
    \psblock(6,0){delay}{$z^{-1}$}
    \pscircleop[oplength=0.25,operation=times](4,0){op3}
    \pssignal(3,0){coefy}{$\frac{15}{16}$}
    \nclist{->}{ncline}{x,op1,op2,y}
    \nclist{->}{ncline}{coefx,op1}
    \nclist{->}{ncline}{ydot,delay,op3,op2}
    \nclist{->}{ncline}{coefy,op3}
  \end{pspicture}
\end{center}

--