Discrete-Time Systems II — Difference Equations & Correlation離散時間系統（二）—— 差分方程與相關離散時間システム（II）—— 差分方程式と相関

3.1 Linear Constant-Coefficient Difference Equations3.1 線性常係數差分方程3.1 線形定係数差分方程式

Chapter 2 described an LTI system by its impulse response $h[n]$ and the convolution sum. A second — and computationally far cheaper — description is the linear constant-coefficient difference equation（線性常係數差分方程, LCCDE）: the present output is computed from a finite number of past outputs, the present input, and past inputs.

第 2 章用脈衝響應 $h[n]$ 與卷積和來描述一個 LTI 系統。第二種描述方式——而且計算上便宜得多——是線性常係數差分方程（linear constant-coefficient difference equation，LCCDE）：把當前的輸出用有限個過去輸出、當前輸入與過去輸入算出來。

第 2 章では、LTI システムをそのインパルス応答 $h[n]$ と畳み込み和によって記述した。もう一つの——そして計算量がはるかに小さい——記述が線形定係数差分方程式（linear constant-coefficient difference equation, LCCDE）である。現在の出力を、有限個の過去の出力と、現在および過去の入力から計算する。

3.1.1 General form3.1.1 一般形式3.1.1 一般形

Solving for the present output gives the recursive computational form:

把當前輸出解出來，就得到遞迴計算形式：

現在の出力について解くと、再帰的な計算形が得られる：

• $N$ = order of the equation（階數）= number of past outputs that feed back.
• Linear: $y$ and $x$ terms appear to the first power only, never multiplied together.
• Constant-coefficient: the $a_k, b_k$ do not depend on $n$ — this is what makes the system time-invariant.

• $N$ ＝ 方程的階數（order）＝ 回授進來的過去輸出個數。
• 線性：$y$ 與 $x$ 的各項都只出現一次方，彼此也不相乘。
• 常係數：$a_k$、$b_k$ 不隨 $n$ 改變——這正是系統非時變的原因。

• $N$ ＝ 方程式の次数 ＝ 帰還してくる過去の出力の個数。
• 線形：$y$ と $x$ の項は一次のみで現れ、互いに掛け合わされることはない。
• 定係数：$a_k, b_k$ は $n$ に依存しない——これが系を時不変にしている。

3.1.2 Recursive vs non-recursive systems3.1.2 遞迴系統與非遞迴系統3.1.2 再帰系と非再帰系

A classic non-recursive example is the $L$-point moving average（移動平均）, used everywhere for smoothing noisy data:

最經典的非遞迴例子是$L$ 點移動平均（moving average），到處都用它來平滑含雜訊的資料：

非再帰系の古典的な例が$L$ 点移動平均（moving average）であり、雑音を含むデータの平滑化にいたるところで使われる：

3.2 Solving LCCDEs3.2 差分方程的解法3.2 差分方程式を解く

The direct method（直接法）solves an LCCDE in the time domain by splitting the total solution into two parts:

直接法（direct method）在時域裡解差分方程，做法是把全解拆成兩部分：

直接法（direct method）は、全解を次の二つに分けて時間領域で差分方程式を解く：

3.2.1 Complementary (homogeneous) solution3.2.1 齊次解（互補解）3.2.1 斉次解（相補解）

The complementary solution（齊次解 / 互補解） $y_c[n]$ solves the homogeneous equation — the LCCDE with the input forced to zero:

齊次解／互補解（complementary solution）$y_c[n]$ 解的是齊次方程——也就是把輸入強迫設為零的那條差分方程：

斉次解（相補解、complementary solution） $y_c[n]$ は斉次方程式——入力を強制的に零とした差分方程式——の解である：

Procedure. Assume an exponential trial solution $y[n] = \lambda^{\,n}$. Substituting and dividing out $\lambda^{\,n-N}$ yields the characteristic polynomial（特徵多項式）:

步驟。先假設指數型試解 $y[n] = \lambda^{\,n}$。代入後約掉 $\lambda^{\,n-N}$，就得到特徵多項式（characteristic polynomial）：

手順。指数型の試行解 $y[n] = \lambda^{\,n}$ を仮定する。代入して $\lambda^{\,n-N}$ で割ると、特性多項式（characteristic polynomial）が得られる：

With $N$ distinct roots $\lambda_1, \dots, \lambda_N$ (the natural modes / poles of the system):

若有 $N$ 個相異根 $\lambda_1, \dots, \lambda_N$（系統的自然模態／極點）：

相異なる $N$ 個の根 $\lambda_1, \dots, \lambda_N$（系の固有モード／極）をもつ場合：

• Repeated root $\lambda_1$ of multiplicity $m$: contributes $(C_1 + C_2 n + \cdots + C_m n^{m-1})\lambda_1^{\,n}$.
• Complex-conjugate pair $\lambda = r e^{\pm j\theta}$: contributes a real damped oscillation $r^{\,n}(A\cos\theta n + B\sin\theta n)$.
• The constants $C_k$ are not fixed yet — they come from the initial conditions, and only after adding the particular solution.

• 重根 $\lambda_1$、重數為 $m$：貢獻 $(C_1 + C_2 n + \cdots + C_m n^{m-1})\lambda_1^{\,n}$。
• 共軛複數對 $\lambda = r e^{\pm j\theta}$：貢獻一個實數的阻尼振盪 $r^{\,n}(A\cos\theta n + B\sin\theta n)$。
• 常數 $C_k$ 此時還沒定下來——它們由初始條件決定，而且必須等加上特解之後才能求。

• 重複度 $m$ の重根 $\lambda_1$：$(C_1 + C_2 n + \cdots + C_m n^{m-1})\lambda_1^{\,n}$ を寄与する。
• 共役複素根の対 $\lambda = r e^{\pm j\theta}$：実数の減衰振動 $r^{\,n}(A\cos\theta n + B\sin\theta n)$ を寄与する。
• 定数 $C_k$ はこの時点ではまだ決まっていない——初期条件から定まり、しかも特解を加えたあとでなければ決められない。

3.2.2 Particular solution3.2.2 特解3.2.2 特解

The particular solution（特解） $y_p[n]$ is any one solution of the full equation with the input present, valid for $n \ge 0$. Assume a form that mimics the input:

特解（particular solution）$y_p[n]$ 是完整方程（輸入存在時）的任何一個解，在 $n \ge 0$ 成立。假設的形式要模仿輸入的樣子：

特解（particular solution） $y_p[n]$ は、入力がある状態での完全な方程式の解のうち任意の一つで、$n \ge 0$ で成り立つ。入力の形を真似た形を仮定する：

Substitute the assumed form into the LCCDE and solve for $K$ by matching coefficients.

把假設的形式代回差分方程，比對係數即可解出 $K$。

仮定した形を差分方程式に代入し、係数を比較して $K$ を求める。

3.2.3 Total solution recipe3.2.3 全解的解題流程3.2.3 全解の手順

3.2.4 Worked example3.2.4 完整範例3.2.4 例題

3.3 Zero-Input and Zero-State Response3.3 零輸入響應與零狀態響應3.3 零入力応答と零状態応答

The same total response can be split a different way, by cause instead of by mathematical form:

同樣一條全響應還可以用另一種方式拆開——不是依數學形式，而是依「成因」來分：

同じ全応答は、数学的な形ではなく「原因」によって、別のやり方で分解することもできる：

• Zero-input response（零輸入響應） $y_{zi}[n]$: set $x[n] = 0$; the stored initial conditions alone drive the output. It is purely a combination of natural modes $\lambda_k^n$.
• Zero-state response（零狀態響應） $y_{zs}[n]$: set $y[-1]=\dots=y[-N]=0$ (system at rest); the input alone drives the output. This is the part equal to $x[n]*h[n]$, and it contains both natural-mode terms and the particular (forced) term.

• 零輸入響應（zero-input response）$y_{zi}[n]$：令 $x[n] = 0$，只由儲存的初始條件驅動輸出。它純粹是自然模態 $\lambda_k^n$ 的組合。
• 零狀態響應（zero-state response）$y_{zs}[n]$：令 $y[-1]=\dots=y[-N]=0$（系統靜止），只由輸入驅動輸出。這一部分才等於 $x[n]*h[n]$，而且它同時含有自然模態項與特解（受迫）項。

• 零入力応答（zero-input response） $y_{zi}[n]$：$x[n] = 0$ とおき、蓄えられた初期条件だけが出力を駆動する。これは固有モード $\lambda_k^n$ の組合せそのものである。
• 零状態応答（zero-state response） $y_{zs}[n]$：$y[-1]=\dots=y[-N]=0$（系は静止）とおき、入力だけが出力を駆動する。$x[n]*h[n]$ に等しいのはこの部分であり、固有モード項と特解（強制）項の両方を含む。

A classic exam point: a system given non-zero initial conditions is not LTI (the stored state destroys both linearity and time invariance), but its zero-state part $y_{zs}[n]$ is. Only when all initial conditions are zero (the system is at rest) does the convolution relation $y_{zs}[n] = x[n]*h[n]$ hold.

考試常考：給了非零初始條件的系統不是 LTI（線性性與時不變性會被初始狀態破壞），但它的零狀態部分 $y_{zs}[n]$ 才是 LTI——唯有在初始條件全為零（系統靜止）時，系統才滿足 $y_{zs}[n] = x[n]*h[n]$ 這個卷積關係。

試験でよく問われる点：非零の初期条件を与えられた系は LTI でない（蓄えられた状態が線形性も時不変性も壊す）が、その零状態部分 $y_{zs}[n]$ は LTI である。初期条件がすべて零（系が静止）のときに限り、畳み込みの関係 $y_{zs}[n] = x[n]*h[n]$ が成り立つ。

3.4 Impulse Response of LCCDEs3.4 線性常係數差分方程的脈衝響應3.4 線形定係数差分方程式のインパルス応答

3.4.1 Impulse response from a difference equation3.4.1 由差分方程求脈衝響應3.4.1 差分方程式からインパルス応答を求める

The impulse response is the zero-state response to $x[n] = \delta[n]$. For a recursive LCCDE, the simplest route is to iterate the recursion with rest initial conditions; for low orders the closed form follows from the natural modes, since for $n > M$ the input terms vanish and $h[n]$ must be a pure combination of modes:

脈衝響應就是對 $x[n] = \delta[n]$ 的零狀態響應。對遞迴型差分方程而言，最簡單的路徑是在靜止初始條件下迭代遞迴式；階數低時，封閉解可以直接由自然模態寫出來，因為 $n > M$ 之後輸入項全部消失，$h[n]$ 必定是模態的純粹組合：

インパルス応答とは、$x[n] = \delta[n]$ に対する零状態応答のことである。再帰型の差分方程式では、静止初期条件のもとで漸化式を反復するのが最も簡単な道筋である。次数が低ければ閉じた形は固有モードから直ちに得られる。$n > M$ では入力項が消え、$h[n]$ はモードの純粋な組合せでなければならないからである：

3.4.2 FIR vs IIR3.4.2 FIR 與 IIR3.4.2 FIR と IIR

The structure of the LCCDE determines the impulse-response length:

差分方程的結構決定了脈衝響應的長度：

差分方程式の構造が、インパルス応答の長さを決める：

3.5 Stability3.5 穩定性3.5 安定性

A system is BIBO stable（有界輸入有界輸出穩定） if every bounded input produces a bounded output: $|x[n]| \le B_x < \infty$ for all $n$ implies $|y[n]| \le B_y < \infty$ for all $n$.

一個系統是 BIBO 穩定（有界輸入有界輸出穩定）的，如果每一個有界輸入都產生有界輸出：對所有 $n$ 有 $|x[n]| \le B_x < \infty$，就能推出對所有 $n$ 有 $|y[n]| \le B_y < \infty$。

ある系が BIBO 安定（有界入力有界出力安定）であるとは、任意の有界入力が有界出力を生むことをいう。すなわち、すべての $n$ について $|x[n]| \le B_x < \infty$ であれば、すべての $n$ について $|y[n]| \le B_y < \infty$ が成り立つ。

For an LTI system, $|y[n]| = \big|\sum_k h[k]x[n-k]\big| \le B_x \sum_k |h[k]|$, which yields the complete characterization:

對 LTI 系統而言，$|y[n]| = \big|\sum_k h[k]x[n-k]\big| \le B_x \sum_k |h[k]|$，由此得到完整的判準：

LTI 系では $|y[n]| = \big|\sum_k h[k]x[n-k]\big| \le B_x \sum_k |h[k]|$ が成り立ち、ここから完全な判定条件が得られる：

The condition is also necessary: if $S = \infty$, the bounded "worst-case" input $x[n] = \operatorname{sgn}\!\big(h[-n]\big)$ makes $y[0] = \sum_k |h[k]|$ blow up.

這個條件同時也是必要的：若 $S = \infty$，取有界的「最壞情況」輸入 $x[n] = \operatorname{sgn}\!\big(h[-n]\big)$，就會讓 $y[0] = \sum_k |h[k]|$ 發散。

この条件は必要でもある：$S = \infty$ ならば、有界な「最悪の場合」の入力 $x[n] = \operatorname{sgn}\!\big(h[-n]\big)$ をとれば $y[0] = \sum_k |h[k]|$ が発散してしまう。

3.6 Block Diagram Realizations3.6 方塊圖實現3.6 ブロック図による実現

An LCCDE maps directly to hardware/software using three primitive elements（三種基本元件）:

一條差分方程可以用三種基本元件直接對映到硬體／軟體：

差分方程式は、三つの基本素子を使ってハードウェア／ソフトウェアへ直接対応づけられる：

• Unit delay（單位延遲） $z^{-1}$: stores one sample; output is the input one sample ago.
• Multiplier（乘法器）: scales a signal by a constant coefficient.
• Adder（加法器）: sums two or more signals.

• 單位延遲（unit delay）$z^{-1}$：儲存一個樣本；輸出是一個樣本之前的輸入。
• 乘法器（multiplier）：把訊號乘上一個常數係數。
• 加法器（adder）：把兩個以上的訊號相加。

• 単位遅延（unit delay） $z^{-1}$：1 標本を蓄える。出力は 1 標本前の入力である。
• 乗算器（multiplier）：信号を定数係数倍する。
• 加算器（adder）：二つ以上の信号を加え合わせる。

Writing the second-order LCCDE in computational form,

把二階差分方程寫成計算形式，

2 次の差分方程式を計算形で書き、

and drawing one branch per term gives the Direct Form I（直接 I 型） realization — a feed-forward tapped delay line on the input plus a feedback tapped delay line on the output:

再對每一項畫出一條分支，就得到直接 I 型（Direct Form I）實現 —— 輸入端一條前饋的抽頭延遲線，加上輸出端一條回授的抽頭延遲線：

各項ごとに枝を 1 本描くと、直接形 I（Direct Form I）による実現が得られる——入力側のフィードフォワード型タップ付き遅延線に、出力側の帰還型タップ付き遅延線を足したものである：

3.7 Correlation3.7 相關3.7 相関

Correlation（相關） measures the degree of similarity between two sequences as a function of relative time shift. It is the standard tool for detecting a known waveform inside noise (radar, sonar, GPS), estimating time delays, and finding hidden periodicity.

相關（correlation）量測兩條序列在不同相對時間平移下的相似程度。它是在雜訊中偵測已知波形（雷達、聲納、GPS）、估測時間延遲，以及找出隱藏週期性的標準工具。

相関（correlation）は、二つの系列の相似度を相対的な時間シフトの関数として測る。雑音に埋もれた既知波形の検出（レーダ、ソナー、GPS）、時間遅延の推定、隠れた周期性の発見に使われる標準的な道具である。

3.7.1 Cross-correlation and its relation to convolution3.7.1 互相關及其與卷積的關係3.7.1 相互相関と畳み込みとの関係

The integer $\ell$ is the lag（延遲）: $r_{xy}[\ell]$ asks "how much does $x$ look like $y$ shifted right by $\ell$?" Reversing the roles only mirrors the lag axis:

整數 $\ell$ 稱為延遲（lag）：$r_{xy}[\ell]$ 問的是「把 $y$ 右移 $\ell$ 之後，$x$ 和它有多像？」把兩者角色對調，只是把延遲軸鏡射過來：

整数 $\ell$ はラグ（遅れ）と呼ばれる：$r_{xy}[\ell]$ は「$y$ を右へ $\ell$ だけずらしたとき、$x$ はそれにどれだけ似ているか」を問う。両者の役割を入れ替えると、ラグ軸が鏡映されるだけである：

Comparing with the convolution sum shows correlation is convolution with one sequence time-reversed — so all convolution machinery applies:

和卷積和對照就會發現，相關就是把其中一條序列時間反轉後的卷積 —— 因此卷積的所有工具都能套用：

畳み込み和と見比べると、相関とは一方の系列を時間反転したうえでの畳み込みそのものであることが分かる——したがって畳み込みの道具立てがすべて使える：

3.7.2 Autocorrelation properties3.7.2 自相關的性質3.7.2 自己相関の性質

The autocorrelation（自相關） is the correlation of a signal with itself, $r_{xx}[\ell] = \sum_n x[n]\,x[n-\ell]$. Its key properties:

自相關（autocorrelation）是訊號與自己的相關，$r_{xx}[\ell] = \sum_n x[n]\,x[n-\ell]$。它的主要性質：

自己相関（autocorrelation）とは、信号とそれ自身との相関 $r_{xx}[\ell] = \sum_n x[n]\,x[n-\ell]$ のことである。主な性質は次のとおり：

and for the cross-correlation the Cauchy–Schwarz bound gives $\big|r_{xy}[\ell]\big| \le \sqrt{r_{xx}[0]\,r_{yy}[0]} = \sqrt{E_x E_y}$. Dividing by the bound yields the normalized correlation（正規化相關係數）, handy because it is scale-free:

至於互相關，Cauchy–Schwarz 不等式給出界限 $\big|r_{xy}[\ell]\big| \le \sqrt{r_{xx}[0]\,r_{yy}[0]} = \sqrt{E_x E_y}$。除以這個界限就得到正規化相關係數（normalized correlation），很好用，因為它與尺度無關：

相互相関については、Cauchy–Schwarz の不等式から $\big|r_{xy}[\ell]\big| \le \sqrt{r_{xx}[0]\,r_{yy}[0]} = \sqrt{E_x E_y}$ という上界が得られる。この上界で割れば正規化相関係数（normalized correlation）になり、尺度に依存しないので扱いやすい：

3.7.3 Autocorrelation of periodic sequences3.7.3 週期序列的自相關3.7.3 周期系列の自己相関

Periodic (power) signals have infinite energy, so the sum is replaced by a time average. For $x[n]$ periodic with period $N$:

週期（功率）訊號的能量無限，所以要把求和改成時間平均。對週期為 $N$ 的 $x[n]$：

周期（電力）信号はエネルギーが無限なので、和を時間平均で置き換える。周期 $N$ をもつ $x[n]$ については：

3.7.4 Correlation in action3.7.4 相關的實際應用3.7.4 相関の実際

The canonical application is radar/sonar delay estimation（雷達延遲估測）: transmit $x[n]$, receive an attenuated, delayed, noisy echo

最典型的應用是雷達／聲納的延遲估測（delay estimation）：發射 $x[n]$，收到一個被衰減、延遲並混入雜訊的回波

典型的な応用がレーダ／ソナーにおける遅延推定（delay estimation）である：$x[n]$ を送信し、減衰し遅延し雑音の乗ったエコーを受信する

Cross-correlating receive with transmit gives $r_{yx}[\ell] = \alpha\, r_{xx}[\ell - D] + r_{wx}[\ell]$. Since $r_{xx}$ peaks at zero lag and the noise term is small, $r_{yx}[\ell]$ peaks at $\ell = D$ — read off the round-trip delay, hence target range. The same idea underlies GPS code acquisition and pattern matching.

把接收訊號與發射訊號做互相關，得到 $r_{yx}[\ell] = \alpha\, r_{xx}[\ell - D] + r_{wx}[\ell]$。由於 $r_{xx}$ 的峰值在零延遲、雜訊項很小，$r_{yx}[\ell]$ 的峰值就落在 $\ell = D$ —— 直接讀出來回傳播的延遲，也就得到目標距離。GPS 的碼捕獲與樣式比對背後也是同一個想法。

受信信号と送信信号の相互相関をとると $r_{yx}[\ell] = \alpha\, r_{xx}[\ell - D] + r_{wx}[\ell]$ が得られる。$r_{xx}$ はラグ 0 でピークをとり、雑音項は小さいので、$r_{yx}[\ell]$ は $\ell = D$ でピークをとる——往復の遅延、すなわち目標までの距離がそのまま読み取れる。GPS のコード捕捉やパターンマッチングも同じ考え方に支えられている。

3.8 Chapter Summary3.8 本章總結3.8 本章のまとめ

Exercises習題演習問題

All four problems of ICE503 DSP Homework #3, fully worked.

ICE503 DSP 作業 #3 的全部四題，附完整解答。

ICE503 DSP 宿題 #3 の全 4 問を、完全に解いたもの。

clc; clear;

% LCCDE: y[n] + 2y[n-1] - 8y[n-2] = 80 x[n]
a = [1, 2, -8];          % output-side coefficients [a0 a1 a2]
b = 80;                  % input-side coefficient
x = ones(1, 5);          % x[n] = u[n] for 0 <= n <= 4
y_m1 = -8;               % y[-1]
y_m2 = -9;               % y[-2]

% (a) closed-form solution from Exercise 2
n = 0:4;
y_closed = 16*(-4).^n + 24*2.^n - 16;
disp('(a) closed form:');  disp(y_closed);

% (b) direct recursion with a for loop
y = zeros(1, 5);
for k = 1:5
    if k == 1          % n = 0: use y[-1], y[-2]
        y(k) = b*x(k) - a(2)*y_m1   - a(3)*y_m2;
    elseif k == 2      % n = 1: use y[0], y[-1]
        y(k) = b*x(k) - a(2)*y(k-1) - a(3)*y_m1;
    else               % n >= 2
        y(k) = b*x(k) - a(2)*y(k-1) - a(3)*y(k-2);
    end
end
disp('(b) for loop:');  disp(y);

% (c) filtic + filter
zi = filtic(b, a, [-8, -9]);     % past outputs [y(-1), y(-2)]
y_filt = filter(b, a, x, zi);
disp('(c) filtic + filter:');  disp(y_filt);

(a) closed form:    24   -32   336   -848   4464
(b) for loop:       24   -32   336   -848   4464
(c) filtic + filter: 24   -32   336   -848   4464

clc; clear;

% (a) generate the windowed cosine on -10 <= n <= 10
n = -10:10;
x = zeros(1, length(n));
for i = 1:length(n)
    if abs(n(i)) <= 5
        x(i) = cos((1/5)*pi*n(i));
    end                          % else: stays 0
end

% (b) autocorrelation r_xx[l] = sum_n x[n] x[n-l]
lags = -(length(n)-1) : (length(n)-1);
r = zeros(1, length(lags));
for k = 1:length(lags)
    s = 0;
    for i = 1:length(x)
        j = i + lags(k);                 % index of x[n - l] (shifted copy)
        if j >= 1 && j <= length(x)      % keep index inside the array
            s = s + x(i) * x(j);
        end
    end
    r(k) = s;
end
% one-line check: r = xcorr(x);  gives the same 41 values

figure;
subplot(2,1,1);
stem(n, x, 'filled'); grid on;
title('x[n] = cos(\pi n/5),  |n| \leq 5'); xlabel('n'); ylabel('x[n]');
subplot(2,1,2);
stem(lags, r, 'filled'); grid on;
title('Autocorrelation r_{xx}[l]'); xlabel('lag l'); ylabel('r_{xx}[l]');