類比轉數位 (A/D Conversion)#

訊號從連續時間訊號轉為類比訊號需要經過兩個步驟，分別是取樣(Sampling)以及量化(Quantization)：

• 取樣（時間離散化）：將 $x(t)$ 轉換為 $x[n]$，決定了訊號的頻率範圍
  • 取樣週期($T$)：兩次取樣之間的時間間隔
  • 取樣率/頻率($f_s$)：$f_s = \frac{1}{T}$
  • 關係式：$t = n \cdot T$
• 量化（數值離散化）：將連續的振幅值轉為有限精度的數位資料，決定了訊號的動態範圍與細節

任意序列的表示 (Representation of Arbitrary Sequence)#

任何離散序列 $x[n]$ 都可以表示為一組加權且延遲的單位脈衝之和 ：

寫成一般式（$a_k$ 即為 $x[k]$ 的值）：

$u[n]$ 與 $\delta[n]$ 的關係

• 累加關係： $u[n]=\sum_{k=0}^\infty \delta[n-k]$ 或 $u[n]=\sum_{k=-\infty}^n \delta[k]$
• 差分關係： $\delta[n]=u[n]-u[n-1]$

複數平面與尤拉公式 (Euler’s Formula)#

在處理複數指數訊號時常用到：

• 尤拉公式：$e^{j \phi}=\cos \phi+j\sin \phi$
• 極座標表示：$A=|A| e^{j \phi_1}$，$\alpha= |\alpha| e^{j \phi_2}$
• 因此 $x[n]=A \alpha^n= |A| |\alpha|^n e^{j(n \phi_2 + \phi_1)}$

✍️TODO 以後想到再補圖

基本定義#

假設一個連續時間訊號 $x(t) = \cos(\omega t)$

• $\omega$：角頻率 (Radian Frequency)，單位：rad/s（徑度/秒）
• $t$：連續時間，單位：sec（秒）

NOTE

複習一下高中內容：在圓周運動或簡諧運動的數學推導中，如果使用「圈數」，公式裡會不斷出現 $2 \pi$ 這個係數，顯得累贅。為了讓方程式更簡潔，我們直接把 $2 \pi$ 乘進去，定義出角頻率（以 $\omega$ 表示）

取樣過程（Sampling Process）#

1. 每隔 $T$ 秒取樣一次，$T$: 取樣週期 (Sampling Period)

   $$\begin{aligned} {\text{時間 } t\text{：}} & 0, T, 2T, 3T, \ldots \\ t &= nT \\ n &= 0, 1, 2, 3, \ldots \end{aligned}$$

2. 離散化

   $$\begin{aligned} x(t) &= \cos(\omega {\color{#3071c4}t}) \\ &\downarrow {\color{#3071c4}t = nT} \\ x({\color{#3071c4}nT}) &= \cos(\omega {\color{#3071c4}nT}) \\ &\Downarrow \\ x[n] &= \cos(\omega nT) \end{aligned}$$

3. 取樣頻率轉換 取樣頻率：$f_s = \frac{1}{T}$ Hz，取樣角頻率為 $\omega_s = \frac{2 \pi}{T}$ rad/s

   TIP

   $\omega = 2 \pi f \xrightarrow{f = \frac{1}{T}} \omega = \frac{2 \pi}{T}$

離散時間的週期性 (Periodicity in Discrete Time)#

假設這裡有兩個訊號：$x$ 和 $x_1$

很明顯這是兩個不一樣頻率的訊號，但如果我們對它們進行採樣，會發現一個很有趣的現象：

可以觀察到雖然 $x$ 跟 $x_1$ 不同，但取樣後都是 $cos(\omega nT)$，$\therefore x[n] = x_1[n]$。

NOTE

這裡再假設一個 $x_2(t) = \cos((-\omega+{\color{#3071c4}\omega_s})t)$，取樣後會得到 $x_2[n] = \cos(-\omega nT)$，但因為 cosine 是偶函數，所以 $x_2[n] = \cos(-\omega nT) = \cos(\omega nT) = x[n]$

我們可以得到一個結論，在離散的世界裡，很多不同頻率的連續波，採樣後居然會看起來一模一樣，這個就是所謂的 「混疊現象」(Aliasing Phenomenon)。

🔗圖片來源：https://www.ni.com/docs/zh-TW/bundle/labwindows-cvi/page/advancedanalysisconcepts/aliasing.html

完美還原訊號的數學條件#

Ex.1#

Determine which input signals to a digital filter or DSP system will be aliased by the given period $\text{Т}$?

1. $x(t) = 2\cos(10t)\text{ , T = 0.1s}$
2. $x(t) = 8\cos(15t)\text{ , T = 0.2s}$

ANS

1. $\omega_s = \frac{2\pi}{T} = \frac{2\pi}{0.1} = 20\pi$
   Nyquist Frequency: $\frac{\omega_s}{2} = 10\pi\approx 31.4$
   $\therefore \omega_N = 10 < \frac{\omega_s}{2}$，沒有混疊
2. 同理，$\frac{\omega_s}{2} = \frac{\frac{2\pi}{0.2}}{2} = 5\pi\approx 15.7$
   $\therefore \omega_N = 15 < \frac{\omega_s}{2}$，沒有混疊

Ex.2#

Determine if the following signals will be aliased. If the signal is aliased into having the same sample value as a lower frequency sinusoidal signal, determine the lower sinusoidal signal.

1. $x(t) = 7\cos(25t)\text{ , T = 0.1s}$
2. $x(t) = 3\cos(160t)\text{ , T = 0.02s}$

ANS

1. $\frac{\omega_s}{2} = \frac{\frac{2\pi}{0.1}}{2} = 10\pi\approx 31.4$
   $\therefore \omega_N = 25 < \frac{\omega_s}{2}$，沒有混疊
2. $\frac{\omega_s}{2} = \frac{\frac{2\pi}{0.02}}{2} = 50\pi\approx 157$
   $\therefore \omega_N = 160 > \frac{\omega_s}{2}$，有混疊
   由於混疊高頻摺疊至低頻的現象，$-\omega_1 + \omega_s = 160$
   又 $\omega_s = 100\pi \Rightarrow \omega_1 = 100\pi - 160$
   $\therefore x_1(t) = 3\cos((100\pi - 160)t) \approx 3\cos(154.2t)$

奈奎斯特頻率 / 摺疊頻率（Nyquist / Folding frequency）#

低通濾波器規格 (Lowpass Digital Filter Specification)#

• 概念圖
• 規格定義圖

  NOTE

  在工程設計中，只要增益曲線（綠線）不要碰到「容差邊限」或稱「禁制區」（黑色的斜線區域）就是合格的濾波器了

高通濾波器規格 (Highpass Digital Filter Specification)#

• 概念圖
• 規格定義圖

帶通濾波器規格 (Bandpass Digital Filter Specification)#

• 概念圖
• 規格定義圖

帶阻濾波器規格 (Bandstop Digital Filter Specification)#

• 概念圖
• 規格定義圖

  TIP

  某些帶阻濾波器又被稱為帶陷濾波器 (Notch Filter)，Notch 是凹口或是刻痕的意思，一般來說，Notch Filter 專指頻帶非常窄的帶阻濾波器。

Ex.1#

Draw the graphical specification for a digital highpass filter out to the folding frequency in rad/s that will not change the gain above 500 rad/s by more than ±3 dB, while reducing the gain below 200 rad/s by more than 40 dB. The sampling time T = 0.001 s.

ANS

1. 先確定 x 軸的盡頭，也就是 $\omega_f$ 的部分，由於取樣時間 $T = 0.001\text{s}$
   $\because\text{Nyquist Frequency} = \frac{\pi}{T}$ 💡為什麼？
   $\therefore \omega_f = \frac{\pi}{0.001} = 1000\pi$
2. 畫出阻帶，題目說 200 rad/s 以下增益要衰減 40 dB，因此阻帶容許最大增益為 -40 dB，即 $g_{smax} = -40\text{ db}$
3. 畫出通帶，題目說 500 rad/s 以上增益改變不能超過 ±3 dB，以無損增益 0 dB 為基準，可得知允許的範圍最高 $g_{pmax} = 3\text{ dB}$，最低 $g_{pmin} = -3\text{ dB}$

規格圖如下所示：

Ex.2#

Draw the graphical specification for a bandstop digital filter out to its folding frequency that will reduce the gain between 1,000 rad/s and 5,000 rad/s by more than 60 dB, but will not reduce the gain above 10,000 rad/s or below 150 rad/s by more than 3 dB. The sampling rate is 10,000 Hz.

ANS

$\text{Nyquist Frequency} = \frac{\omega_s}{2} = \frac{\cancel{2}\pi\cdot f_s}{\cancel{2}} = 10000\pi$
後面跟上一題類似故省略
規格圖如下所示：

正規化角頻率 (Normalized Angular Frequency)#

在處理離散訊號時，為了避免公式裡一直帶著取樣週期 $T$ 而顯得臃腫，我們通常會定義一個新的變數：

這裡的 $\Omega = \omega T$ 稱為 離散角頻率 或 正規化角頻率，單位為 rad。

1. 單位樣本/脈衝訊號 (Unit Impulse)#

根據定義，$\delta[n]$ 只在 $n=0$ 時有值，其餘皆為 0。

2. 單位步階訊號 (Unit-step Signal)#

對於 $u[n]$，訊號在 $n \ge 0$ 時恆為 1。

REMARK

當 $|r| < 1$ 時，無窮等比級數之和為： $1 + r + r^2 + r^3 + \dots = \frac{1}{1-r}$

將 $u[n]$ 代入 Z 轉換公式，令 $r = z^{-1}$：

Z 轉換的平移特性#

這裡先提一下 Z 轉換的平移特性，在時域延遲 $k$ 個單位的 $x[n-k]$，經過 Z 轉換後會乘上 $z^{-k}$，即為 $z^{-k}\cdot X(z)$

舉個栗子|⩊･)ﾉ🌰#

有個 DSP 系統，它的差分方程式為：$y[n] -2\cdot y[n-1] = 0.5\cdot x[n-1]$，求它的 Block Diagram 。
畫出來要像下面這樣：

我們對差分方程式進行整理：

於是我們就能畫出：

轉移函數#

系統的轉移函數定義為輸出 $Y(z)$ 與輸入 $X(z)$ 的比值。承接上面的例子，首先對差分方程式進行 Z 轉換，根據上面證明的平移特性，我們可以很快的寫出：

這個轉移函數 $T(z)$ 是描述系統特性的核心，它決定了系統面對任何輸入時的反應：

只要有了 $T(z)$，我們就能輕鬆算出任何輸入對應的結果：$Y(z) = T(z)\cdot X(z)$

振幅響應與增益 (Magnitude Response / Gain)#

• 系統的增益即為轉移函數在 $z = e^{j\Omega}$ 時的絕對值大小：$\text{Gain} = |T(z)|_{z=e^{j\Omega}}$
• 若要以分貝表示，公式為：$\text{Gain}_{\text{dB}} = 20\cdot\log(|T(z)|_{z=e^{j\Omega}})$
• 可利用尤拉公式 (Euler formula) $e^{j\Omega} = \cos(\Omega) + j\sin(\Omega)$ 將數值代入轉移函數，計算特定頻率下的增益

Remark: Euler Formula

$$\begin{aligned} e^{j\theta} &= \cos(\theta) + j\sin(\theta), \quad j = \sqrt{-1} \\ \therefore z = e^{j\Omega} &= \cos(\Omega) + j\sin(\Omega) \end{aligned}$$

代入轉移函數後會得到一個複數 $T(e^{j\Omega}) = a + jb$，其絕對值（也就是增益）算法為： $|T(e^{j\Omega})| = \sqrt{a^2 + b^2}$

相位響應 (Phase Response)#

除了振幅（增益）會改變之外，訊號經過系統後，波形的「相位」也會產生偏移。 將 $z = e^{j\Omega}$ 代入轉移函數後得到的複數 $a + jb$，我們也可以求出它的角度，這個角度就是相位的變化量 $\phi$：

極點與零點 (Poles & Zeros)#

轉移函數 $T(z)$ 通常可以表示成一個帶有分子跟分母的有理函數 (Rational function)： $T(z) = \frac{Y(z)}{X(z)} = \frac{N(z)}{D(z)}$

• 零點 (Zeros)： 讓分子 $N(z) = 0$ 的 $z$ 值。在這些特定的 $z$ 值下，系統的輸出為 0（增益為 0），代表這些頻率的訊號會被完全阻擋。
• 極點 (Poles)： 讓分母 $D(z) = 0$ 的 $z$ 值。在這些特定的 $z$ 值下，系統增益會趨近於無限大，極點的位置對於判斷一個系統是否「穩定」非常關鍵。

轉移函數 (Transfer Function)#

在時域看這長長一串方程式有點痛苦，所以我們把它丟進 Z 轉換的魔法陣裡。利用上一章提到的平移特性（$x[n-k] \xrightarrow{\mathcal{Z}} z^{-k}X(z)$），我們可以把差分方程式整理成系統的轉移函數 $H(z)$：

經過移項與提出公因數後，我們就能得到系統的轉移函數 $H(z)$：

NOTE

這裡的 $H(z)$ 跟前面提到的 $T(z)$ 基本上是同一個東西，其實在 DSP 中用 $H(z)$ 應該是比較普遍的，之後不管看到 $T(z)$ 還是 $H(z)$ 都要知道是在說轉移函數（系統函數）。

舉個栗子|⩊･)ﾉ🌰#

假設我們已知一個一階類比低通濾波器（直流增益為 1，截止頻率 $\omega_c = 10\text{ rad/s}$），其轉移函數為 $T(s) = \frac{10}{s+10}$。 現在我們想要使用脈衝不變法來找它的數位濾波器差分方程式，假設取樣週期 $T = 0.1\text{s}$。

Step 1. 求出連續時間脈衝響應
對 $T(s)$ 取反拉普拉斯轉換：

Step 2. 代入離散時間進行取樣
將 $t = nT = 0.1n$ 代入上述方程式：

Step 3. 進行 Z 轉換
將取樣後的 $y[n]$ 進行 Z 轉換求出 $Y(z)$：

Step 4. 進行增益補償並求出 $H(z)$
為了補償取樣過程中產生的 $1/T$ 倍頻譜放大效應，須將 Z 轉換結果乘上取樣週期 $T$ 以校正增益，使其與原類比濾波器量級一致：

為什麼要乘上取樣週期 $T$？

如果不做增益補償，系統會怎麼樣？我們實際拿上面這個例子算一下「直流增益」（也就是頻率為 0 時的倍率，代入 $s=0$ 或 $z=1$）：

• 原本的類比濾波器： $T(s) = \frac{10}{0+10} = 1$ （訊號 1 倍無損通過）
• 未乘 T 補償的數位濾波器： $Y(z) = \frac{10}{1-e^{-1}} \approx 15.8$ （訊號直接被暴增近 16 倍，系統大爆走）
• 乘上 $T\ (0.1)$ 補償後： $H(z) = \frac{10}{1-e^{-1}} \cdot 0.1 \approx 1.58$

雖然脈衝不變法因為天生的「高頻混疊」缺陷，導致補償後無法完美回到 1，但至少乘上 $T$ 後，已經把增益量級拉回正確的基準線了！

Step 5. 轉回差分方程式
有了轉移函數，我們將其還原為輸入與輸出的關係：

最後反 Z 轉換回時域，大功告成！這就是我們要的 IIR 濾波器：

✍️TODO：圖之後補

進階範例|⩊･)ﾉ🥥：二階系統的脈衝不變法#

Find the digital IIR filter for the second-order lowpass analog filter given by the following equation for $T(s) = \frac{15}{s(s+15)}$ for a sampling period $T = 0.05 \text{ s}$.

Step 1. 部分分式展開 (Partial Fraction Expansion, PFE)
為了方便進行反拉普拉斯轉換，我們先把 $T(s)$ 拆開：

經過計算可得 $A=1, B=-1$。
所以 $T(s)={\color{#0D9488}\frac{1}{s}}-{\color{#D97706}\frac{1}{s+15}}$。

Step 2. 轉換至時域
查閱拉普拉斯轉換表，我們可以得到連續時間的脈衝響應：$y(t) = {\color{#0D9488}1} - {\color{#D97706}e^{-15t}}$

Step 3. 進行取樣與 Z 轉換
代入 $t=nT$，並進行 Z 轉換求出 $T(z)$（原理與一階系統相同，這裡快速帶過）。最終整理後可得到系統的轉移函數：

Step 4. 轉回差分方程式
將轉移函數 $T(z) = \frac{Y(z)}{X(z)}$ 展開並交叉相乘：

接著利用 Z 轉換的平移特性，將等式反轉換回時域，就能得到最終的差分方程式：

舉個栗子|⩊･)ﾉ🌰#

Finding the Step Invariant IIR filter to approximate a first-order lowpass Butterworth (analog) filter $T(s) = \frac{10}{s+10}$, assuming sampling time $T = 0.1 \text{ s}$.

Step 1：求類比系統的步階響應 $Y(s)$
為了求步階響應，我們假設輸入 $X(s)$ 為單位步階函數 $\frac{1}{s}$