Kalman Filter (1)

(這篇文章拖了有點久，JK同學已經在催了....) 想要瞭解Kalman Filter最好還是應該先瞭解線性系統的基本概念，以及一些隨機程序的觀念。這裡先介紹一點很簡單的線性系統理論，特別是狀態空間表示法(State-Space Representation)。 一個離散時間線性系統可以狀態空間方程式表示如下： \[ \begin{array}{rcl} {\bf x}(k+1) & = & {\bf F}(k){\bf x}(k)+{\bf G}(k){\bf u}(k) \\ {\bf y}(k) & = & {\bf H}(k){\bf x}(k) \end{array} \] 這裡我們討論離散時間的系統，因為Kalman Filter利用離散時間動態系統來解釋較容易理解。上述方程式中，第一個式子稱為狀態方程式(state equation)，第二個式子稱為量測方程式(measurement equation)。${\bf y}(k)$為量測輸出，${\bf u}(k)$為系統輸入。對於一個實際的物理系統，我們量測到的輸出(亦即${\bf y}$)，可能只是狀態的一部分，例如一個慣性元件(加速規或者陀螺儀)，我們可以取得的輸出量可能只有載具的位置，角速度、加速度與速度較不容易量測，而實際上若是可以量測得到載具的速度、加速度以及角速度，我們可以對系統做更好的控制。此時，可以利用Kalman所提出來的觀測器的概念來解決這個問題。從數學上來看，載具的加速度、速度以及外力之間的關係，可以寫成一個線性的二階微分方程式：位置(position)的微分是速度、速度的微分是加速度；同樣地，在電路學中，一個RLC電路的電流或者電壓的關係，也可以寫成一個二階微分方程式：以電容儲存的電荷(charge)為例，電荷對時間的微分是電流，而電流對時間的微分又正比於電感上的電壓。從數學上的關係來看，若是電荷量可以透過適當的比例因子與載具的位置建立一個正比的關係，例如：position = k $\times$ charge，k為常數；利用常係數微分方程式的線性關係，可知載具的速度與加速度，和RLC電路中的電流以及電壓也會存在正比的關係。換句話說，我們可以從RLC電路中，得知載具的速度以及加速度的資訊(或者至少可以從"數學上"得知這些訊息)。很明顯的，在電路中量測電荷量、電流以及電壓要簡單多了，因此可以將這個RLC電路，看成載具中慣性元件的一種"觀測器"(observer or estimator)。了解觀測器的概念之後，我們接著來看在數學上如何設計觀測器。

圖 1. 觀測器架構 

從上面的說明，不難發現，在觀測器中有一部分的動態方程式，必須是從"受控體"(plant)中複製過來的；另外，受控體的輸入以及量測輸出的資訊，也必須同時提供給觀測器；因此，不妨假設觀測器的狀態空間方程式如下(參考圖1)：
\[ \begin{array}{rcl} \hat{\bf x}(k+1) & = & {\bf F}(k)\hat{\bf x}(k)+{\bf G}(k){\bf u}(k) + {\bf L}(k)({\bf y}(k)-\hat{\bf y}(k)) \\ \hat{\bf y}(k) & = & {\bf H}(k)\hat{\bf x}(k) \end{array} \] 其中$\hat{\bf x}$代表觀測器的狀態，用以估測實際的系統狀態${\bf x}$，${\bf L}$則是一個待決定的參數，稱為觀測器增益向量(observer gain vector)。定義觀測器誤差為：$\delta{\bf x}(k) := {\bf x} - \hat{\bf x}$，將前面的兩組方程式相減之後，可得到觀測誤差的動態方程式為：
\[ \delta{\bf x}(k+1) = ({\bf F}(k) - {\bf L}(k){\bf H}(k))\delta{\bf x}(k). \] 從前面介紹的觀測器的概念來看，我們當然希望觀測誤差為零，或者至少在經過足夠的時間之後，誤差為零；從上面這個動態方程式來看，$\delta{\bf x}$最後會趨近於零的充分且必要條件為：${\bf F}(k) - {\bf L}(k){\bf H}(k)$這個矩陣必須是漸近穩定的(asymptotically stable)。因此，觀測器的問題現在變成是要找到觀測器的增益向量${\bf L}$使得${\bf F}(k) - {\bf L}(k){\bf H}(k)$為漸近穩定。Kalman告訴我們，是否找的到這樣的向量，與系統的可觀測性(Observability)有關。這邊，我們暫時不討論可觀測性的問題，有興趣的讀者可以參考線性系統的書籍，大部分都有介紹。假設找的到這樣的矩陣${\bf L}$，則可確保由觀測器所得到之狀態空間估測值$\hat{\bf x}$，最後會"追"上系統真正的狀態空間變數${\bf x}$。

但是，若系統有雜訊的話，會發生甚麼問題嗎?由於雜訊多半僅能以機率的方式來描述，從數學上來看，觀測器的狀態無法真正追上受控系統的狀態，狀態空間變數的誤差$\delta{\bf x}(k)$不會趨近於零，而且因為系統有雜訊，$\delta{\bf x}(k)$會是一個隨機程序。我們知道，對一個隨機程序，最好的描述方式是用統計平均值，包括mean, variance等。既然現在無法令$\delta{\bf x}(k)$趨近於零，那在"最佳"的情況下，我們能做到甚麼程度呢？從統計的觀點來看，一個直覺地評估我們設計出來的觀測器的好壞的方式，第一個當然希望$\delta{\bf x}(k)$期望值為零(也就是說，至少可以"期望"${\bf x}$與$\hat{\bf x}$長時間下來的平均值是一樣的)，第二個可能就是希望$\delta{\bf x}(k)$每次估測出來，"跳動"不會太大，在統計學的術語上就是"標準差"(standard deviation，其平方稱為變異量，variance)。這樣應該可以稱為是"最佳"的觀測器。選擇觀測器增益向量${\bf L}$，使得估測誤差$\delta{\bf x}$的mean為零，variance為最小，就是所謂的Kalman Filter。

(待續)

通訊原理 -- 第一章 訊號與系統

訊號與系統

GPS接收機原理上課講義 (I)

Lecture 01

Two-Port Networks

雙埠網路(Two-Port Network)在電子學中佔有重要的地位，大學的基礎電子學課程一直在重複的應用雙埠網路的理論。在電路中，電流可以從端點(terminal)流進或者流出，成對的端點即可構成一個埠(port)。電路中的三個基本元件：電阻、電感以及電容都是雙端點元件，可以形成一個單埠網路(one-port network)，大學的基礎電路學大部份的時間都是在研究單埠網路。電阻、電容以及電感這三個基本元件，稱為被動元件(passive elements)。

圖1. (a)單埠網路 (b)雙埠網路

在電子學中，除了雙端元件之外，還多了三端(例如電晶體)以及四端(例如運算放大器、MOS)元件。此時，電路的分析便需要利用雙埠網路。一般而言，一個複雜的電路(或者網路)可以有n個埠，但是在初學的階段，雙埠網路已經足夠。如圖1(b)所示，雙埠網路有兩組端點對(terminal pairs)，對網路的存取都是透過端點來進行，電流從成對端點的一頭流入，從另一頭流出(當然必須滿足電路的特性，亦即流進的電流必須等於流出的電流)。雙埠網路的目的實際上是希望簡化電路分析的過程，將複雜的網路化簡成一個線性電路，僅需要幾個參數即可描述電路的特性。另外，若要將雙埠網路加入至另一個更大更複雜的網路，此時可將雙埠網路視為一個"黑盒子"(black box)，對於外部更大的網路而言，不需要知道其內部實際電路為何，僅需要知道雙埠輸入輸出之電壓電流關係即可進行電路分析。

雙埠網路的另一個重要特性是"線性化"，也就是將一個複雜的電路，利用簡單的線性模型來表示。在圖1(b)中，$(V_1, I_1)$與$(V_2, I_2)$之間可能存在非常複雜的非線性關係(例如在半導體元件中，電流電壓的關係通常需要高度非線性的方程式來描述)，分析起來不容易。在一些合理的假設條件之下(學過電子學的同學應該知道，就是所謂的小信號條件，small-signal condition)，可以使用線性模型來描述$V_1,V_2,I_1,I_2$之間的關係。線性模型實際上就是一個線性映射(linear transformation)，從$V_1,V_2,I_1,I_2$之中選擇兩個量為自變數，另外兩個量為應變數，則描述這四個量之間數學關係的線性映射為一個$2\times 2$的矩陣，換句話說需要四個參數即可完整描述雙埠網路的線性模型。依據自變數與應變數選擇方式的不同，可得到四種雙埠網路參數，分別敘述如下。

圖2. 電壓源驅動雙埠網路

Impedance Parameters
考慮如圖2所示之雙埠網路，我們稱左邊的埠為Port 1(or input port)，右邊的埠為Port 2 (or output port)。在假設該電路為線性網路的情況下，埠1及埠2的端點電壓與電流之間的關係可以表示為：
\[\begin{array}{l}

{V_1} = {z_{11}}{I_1} + {z_{12}}{I_2}\\
{V_2} = {z_{21}}{I_1} + {z_{22}}{I_2}
\end{array}\hspace{2cm}\mbox{(1)}\]
或者以矩陣的符號可以表示為：

\[\left[ {\begin{array}{*{20}{c}}
{{V_1}}\\
{{V_2}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{z_{11}}}&{{z_{12}}}\\
{{z_{21}}}&{{z_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{I_1}}\\
{{I_2}}
\end{array}} \right] = {\bf{Z}}\left[ {\begin{array}{*{20}{c}}
{{I_1}}\\
{{I_2}}
\end{array}} \right]\]
其中$\{z_{11},z_{12},z_{21},z_{22}\}$稱為雙埠網路的阻抗參數(impedance parameters)，或者稱為$z$參數。 那麼這些參數如何求得？從數學的角度來看，當然是利用偏微分：
\[{z_{11}} = \frac{{\partial {V_1}}}{{\partial {I_1}}},{z_{12}} = \frac{{\partial {V_1}}}{{\partial {I_2}}},{z_{21}} = \frac{{\partial {V_2}}}{{\partial {I_1}}},{z_{22}} = \frac{{\partial {V_2}}}{{\partial {I_2}}}.\]
可是，這樣的物理意義並不明顯。另一種方式則是分別令$I_1=0$(input port open-circuited)或者$I_2=0$(out port open-circuited)，則由(1)式可得
\[{z_{11}} = {\left. {\frac{{{V_1}}}{{{I_1}}}} \right|_{{I_2} = 0}},{z_{12}} = {\left. {\frac{{{V_1}}}{{{I_2}}}} \right|_{{I_1} = 0}},{z_{21}} = {\left. {\frac{{{V_2}}}{{{I_1}}}} \right|_{{I_2} = 0}},{z_{22}} = {\left. {\frac{{{V_2}}}{{{I_2}}}} \right|_{{I_1} = 0}}.\hspace{1cm} \mbox{(2)}\]
這樣的表示方式，應該可以明顯看出$z$參數的物理意義：由於電路中，電流為零代表開路，因此上式的意思係指我們只要適當的將port 1或者port 2開路，並且量測得到對應的電壓電流的關係即可求得$z$參數。

• $z_{11}$ = open-circuit input impedance
• $z_{12}$ = open-circuit transfer impedance from port 1 to port 2
• $z_{21}$ = open-circuit transfer impedance from port 2 to port 1
• $z_{22}$ = open-circuit output impedance

圖3. 決定$z$參數：(a)求$z_{11}$與$z_{21}$, (b)求$z_{21}$與$z_{22}$.

由於$z$參數係將輸出或者輸入開路所得到的，因此也稱為開路阻抗參數(open-circuit impedance parameters)。依據第(2)式，將電壓$V_1$接至port 1，且將port 2開路，並測得$I_1$以及$V_2$(例如利用電流計與電壓計)，如圖3(a)，則可求得$z_{11}$與$z_{21}$為
\[{z_{11}} = \frac{{{V_1}}}{{{I_1}}},{z_{21}} = \frac{{{V_2}}}{{{I_1}}}.\]
同理，欲求得$z_{21}$與$z_{22}$則將電壓$V_2$接至port 2，且將port 1開路，如圖3(b)，並測得$I_2$以及$V_1$，則
\[{z_{12}} = \frac{{{V_1}}}{{{I_2}}},{z_{22}} = \frac{{{V_2}}}{{{I_2}}}.\]

Admittance Parameters
第二種模型，係將兩端的電流，以兩端的電壓來表示，因此稱為導納參數(admittance parameters)。方程式如下：
\[\begin{array}{l}
{I_1} = {y_{11}}{V_1} + {y_{12}}{V_2}\\
{I_2} = {y_{21}}{V_1} + {y_{22}}{V_2}
\end{array}\hspace{2cm}\mbox{(3)}\]
或者以矩陣的符號可以表示為：

\[\left[ {\begin{array}{*{20}{c}}
{{I_1}}\\
{{I_2}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{y_{11}}}&{{y_{12}}}\\
{{y_{21}}}&{{y_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{V_1}}\\
{{V_2}}
\end{array}} \right] = {\bf{Y}}\left[ {\begin{array}{*{20}{c}}
{{V_1}}\\
{{V_2}}
\end{array}} \right]\]
其中$\{y_{11},y_{12},y_{21},y_{22}\}$稱為雙埠網路的導納參數，或者稱為$y$參數。如同在計算$z$參數的方式，適當地令$V_1 = 0$或者$V_2 = 0$，可以求得這些參數，如下：

\[{y_{11}} = {\left. {\frac{{{I_1}}}{{{V_1}}}} \right|_{{V_2} = 0}},{y_{12}} = {\left. {\frac{{{I_1}}}{{{V_2}}}} \right|_{{V_1} = 0}},{y_{21}} = {\left. {\frac{{{I_2}}}{{{V_1}}}} \right|_{{V_2} = 0}},{y_{22}} = {\left. {\frac{{{I_2}}}{{{V_2}}}} \right|_{{V_1} = 0}}.\hspace{1cm} \mbox{(4)}\]

圖4. 決定$y$參數：(a)求$y_{11}$與$y_{21}$, (b)求$y_{21}$與$y_{22}$.

電路中，若令電壓為零，表示係將電路短路，因此$y$參數是將電路的輸出或者輸入短路所得到的，也稱為短路導納參數。

• $y_{11}$ = short-circuit input admittance
• $y_{12}$ = short-circuit transfer admittance from port 2 to port 1
• $y_{21}$ = short-circuit transfer admittance from port 1 to port 1
• $y_{22}$ = short-circuit output admittance

Hybrid Parameters
有些時候電路的$z$與$y$參數可能不存在，因此需要另一種參數的模型。第三種參數的模型，是令$V_1$以及$I_2$為相依變數，$I_1$以及$V_2$視為獨立變數，可得到下列方程式：
\[\begin{array}{l}
{V_1} = {h_{11}}{I_1} + {h_{12}}{V_2}\\
{I_2} = {h_{21}}{I_1} + {h_{22}}{V_2}
\end{array}\hspace{2cm}\mbox{(5)}\]
或者以矩陣的符號可以表示為：

\[\left[ {\begin{array}{*{20}{c}}
{{V_1}}\\
{{I_2}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{h_{11}}}&{{h_{12}}}\\
{{h_{21}}}&{{h_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{I_1}}\\
{{V_2}}
\end{array}} \right] = {\bf{H}}\left[ {\begin{array}{*{20}{c}}
{{I_1}}\\
{{V_2}}
\end{array}} \right]\]

此種參數模型稱為混合參數(hybrid parameters)或者簡稱為$h$參數，經常用於描述電晶體等電子電路元件，理想的變壓器也可以利用$h$參數模型來描述。$h$參數的計算方式如下：

\[{h_{11}} = {\left. {\frac{{{V_1}}}{{{I_1}}}} \right|_{{V_2} = 0}},{h_{12}} = {\left. {\frac{{{V_1}}}{{{V_2}}}} \right|_{{I_1} = 0}},{h_{21}} = {\left. {\frac{{{I_2}}}{{{I_1}}}} \right|_{{V_2} = 0}},{h_{22}} = {\left. {\frac{{{I_2}}}{{{V_2}}}} \right|_{{I_1} = 0}}.\hspace{1cm} \mbox{(6)}\]

圖5. 雙埠網路的$h$參數模型

四個$h$參數的名稱定義如下：

• $h_{11}$ = short-circuit input impedance
• $h_{12}$ = open-circuit reverse voltage gain
• $h_{21}$ = short-circuit forward current gain
• $h_{22}$ = open-circuit output admittance

Inverse Hybrid Parameters
另一組與$h$參數類似的模型，稱為$g$參數，或者反向混合參數 (inverse hybrid parameters)，係將$V_2$以及$I_1$當作相依變數，$V_1$以及$I_2$視為獨立變數，方程式如下：
\[\begin{array}{l}
{I_1} = {g_{11}}{V_1} + {g_{12}}{I_2}\\
{V_2} = {g_{21}}{V_1} + {g_{22}}{I_2}
\end{array}\hspace{2cm}\mbox{(7)}\]
或者以矩陣的符號可以表示為：

\[\left[ {\begin{array}{*{20}{c}}
{{I_1}}\\
{{V_2}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{g_{11}}}&{{g_{12}}}\\
{{g_{21}}}&{{g_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{V_1}}\\
{{I_2}}
\end{array}} \right] = {\bf{G}}\left[ {\begin{array}{*{20}{c}}
{{V_1}}\\
{{I_2}}
\end{array}} \right]\]

圖6. 雙埠網路的$g$參數模型

$g$參數的計算方式如下：

\[{g_{11}} = {\left. {\frac{{{I_1}}}{{{V_1}}}} \right|_{{I_2} = 0}},{g_{12}} = {\left. {\frac{{{I_1}}}{{{I_2}}}} \right|_{{V_1} = 0}},{g_{21}} = {\left. {\frac{{{V_2}}}{{{V_1}}}} \right|_{{I_2} = 0}},{g_{22}} = {\left. {\frac{{{V_2}}}{{{I_2}}}} \right|_{{V_1} = 0}}.\hspace{1cm} \mbox{(8)}\]

四個$h$參數的名稱定義如下：

• $g_{11}$ = open-circuit input admittance
• $g_{12}$ = short-circuit reverse current gain
• $g_{21}$ = open-circuit forward voltage gain
• $g_{22}$ = short-circuit output impedance

其他還有很多種不同的參數模型，例如在高頻微波電路中經常使用的$S$參數模型，就等以後有空再說吧，先這樣。

參考資料：C. K. Alexander and M. N. O. Sadiku, Fundamentals of Electric Circuits, Third Edition.

第7章 GPS載波相位觀測量定位

在第四章中，我們曾經分析電碼相位以及載波相位定位的精確度，在不考慮接收機雜訊以及多路徑效應影響的情況下，電碼相位的定位精確度大多處於公尺等級，而載波相位的定位精確度可達0.01至0.05個週期 (大約2 mm至1 cm)，因此公分等級的精確定位，需要利用載波相位觀測量。載波相位觀測量大多係採用相對定位 (relative positioning) 的方式，其基本概念是將在兩個不同觀測點上的，同一時間的觀測量相減 (亦即差分, difference)，以消去兩者之間的共同誤差，再重新定義問題的參數，以使用這些差分過後的量測量估測這兩個點之間的相對位置向量。GPS系統原本的設計是使用電碼相位觀測量進行定位，GPS系統的設計者並未預期可以提供公分等級的定位精確度，最早展示GPS確實具備精密定位能力，係在1970年代末期，由Counselman以及他在MIT的同事，在論文[1][2]中提出。近代GPS載波相位一次差分以及二次差分的概念，最早都是由Counselman提出。

第7.1節 載波相位與整數未定值

考慮如圖7‑1所示之兩台GPS接收機架構，其天線位置分別以A與B表示，這兩架接收機皆追蹤來自同一顆衛星的載波相位，我們的目的係要精確地量測兩個天線之間的距離d，其中我們假設角度θ₀為已知。依據波動理論，來自衛星的平面波波前 (plane wav front) 代表固定的載波相位，從圖中可以看出來，假設在t₀時某一波前首先到達天線B，則再經過數個整數週期以及一個分數週期 (fractional cycle) 之後，該波前到達天線A，假設該分數週期為Δ₀，由此可知在這兩個天線之間的載波相位觀測量的差必為一特定個數的整數週期加上該分數週期Δ₀。由於正弦波週期信號本身的特性，一開始我們無法判定整數周期的個數，僅能得知該分數周期Δ₀的大小。若以φ_A與φ_B分別表示天線A以及天線B的載波相位，則這兩個天線之間的相位差可以表示為：

\[{\phi _{AB}}({t_0}) = {\phi _A}({t_0}) - \phi ({t_0})=
{\Delta _0} + N\]
其中$N$為未知整數，一般稱為整數週波未定值(integer cycle ambiguity)。若可以決定$N$之值，則由圖7-1(a)中的幾何關係，可以精確地求得$d$之值，如下

(7.1-1)          $d\cos {\theta _0} = \lambda ({\Delta _0} + N).$

若兩個地面接收機的相對位置保持固定，且接收的衛星也固定，則觀測到的載波相位差(Δ₀) 也會維持固定，在此情況下將無從判定N之值。若在一段時間之後，衛星的位置明顯改變，衛星與接收機之間的幾何結構也改變，此時接收機之間的載波相位差也會改變，如圖7‑1(b)所示，此時可利用兩個不同時間所量測到之載波相位差，求解整數週波未定值。假設兩台接收機接連續地追蹤同一顆衛星的載波相位，在時間t₁時，衛星的仰角從θ₀變化至θ₁，兩台接受機所觀測到之相位差也從Δ₀變成Δ₁，若兩台接收機追蹤衛星之信號未曾間斷，則整數週波未定值也不會改變。由於兩台接收機的位置並未改變，因此

(7.1-2)          $d\cos {\theta _1} = \lambda ({\Delta _1} + N).$

將方程式(7.1‑1)與(7.1‑2)重新整理如下

(7.1-3)          $\begin{array}{l} d'\cos {\theta _0} - N = {\Delta _0},\\ d'\sin {\theta _1} - N = {\Delta _1}, \end{array}$

其中d’ = d/λ。注意在上式中單位為周 (cycle)。(7.1‑3)式為二元一次方程組，未知變數為d’與N，其解可以表示為

(7.1-4)          $\left[ {\begin{array}{*{20}{c}}
{d'}\\
N
\end{array}} \right] = \frac{1}{{(\cos {\theta _1} - \cos {\theta _0})}}\left[ {\begin{array}{*{20}{c}}
{ - 1}&1\\
{ - \cos {\theta _1}}&{\cos {\theta _0}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{\Delta _0}}\\
{{\Delta _1}}
\end{array}} \right].$

從上式中不難看出，當θ₀與θ₁很接近時，方程組(7.1‑3)的解，會變成所謂的變態解 (ill-conditioned)，亦即Δ₀與Δ₁產生的微小誤差，會造成方程式解(d’, N)極大的誤差。因此，利用(7.1‑4)式求解整數週波未定值，t₀與t₁的時間間隔必須夠大，以確保θ₀與θ₁差異夠大，此時即可利用(7.1‑4)求解整數週波未定值N以及接收機之間的相對距離d’。這種接收機與衛星之間幾何結構的改變，文獻上稱為幾何多樣性(geometric diversity)。

上述方法雖然在理論上可行，但是為了確保衛星與接收機之間的幾何多樣性，通常需要等待一段時間(大約30分鐘以上)，這在實用上有些時候不是非常方便，因此在1985年時，學者Remondi提出一種可以快速改變衛星-接收機幾何的方法，稱為天線互換 (antenna swap)[4]。假設A與B兩點之間的距離足夠近，在量測到初始相位差Δ₀之後，我們可以將A與B兩個天線位置互換 (注意此時A與B兩台接收機必須持續追蹤同一顆衛星的信號)，若互換時間很短，我們可以假設衛星位置並沒有明顯的改變，則因為天線互換，衛星信號到達天線B的路徑長度會增加d’cosθ₀，反之到達天線A的路徑長度會減少d’cosθ₀。假設在天線交換之後所量測到之相位差為，則交換前後的相位變化與d’cosθ₀之關係可以寫為

\[{\Delta '_0} - {\Delta _0} =  - 2d'\cos {\theta _0}.\]

直接求解可得

\[d' = \frac{{({{\Delta '}_0} - {\Delta _0})}}{{2\cos {\theta _0}}}.\]

代入(7.1‑1)式中，即可求得整數週波未定值為

\[N = d'\cos {\theta _0} - {\Delta _0} =  - \left( {\frac{{{{\Delta '}_0} + {\Delta _0}}}{2}} \right).\]

只要衛星持續追蹤該顆衛星的信號，整數週波未定值也不會改變，換句話說N之值僅需求解一次，一旦求得該衛星的整數未定值，接收機即可自由移動。但是，若接收機並未連續追蹤該顆衛星信號，例如發生週波跳動 (cycle slip) 的現象，則必須重新求解整數週波未定值。

第7.2節 載波相位觀測量與精密定位

在第四章中我們曾介紹GPS的載波相位觀測量，為了方便起見，將載波相位方程式重複如下：

(7.2-1)          $\phi  = {\lambda ^{ - 1}}(r - I + T) + f \cdot (\delta {t_u} - \delta {t^s}) + N + {\varepsilon _\phi },$  (單位為周, cycles)

其中，λ與f分別為載波的波長與頻率，r為衛星與接收機之間的幾何距離，I與T分別為電離層超前 (ionospheric advance) 以及對流層延遲 (tropospheric delay)，衛星以及接收機時鐘偏差分別以δt^s以及δt_u表示，N代表整數週波未定值，表示載波相位量測量的模型誤差。載波相位觀測量與電碼相位觀測輛主要差別有兩點：首先，電碼相位觀測量基本上並沒有整數未定值的問題，而載波相位觀測量因為具有整數未定值，必須先解出未定值之後 (或者消去未定值)，才有辦法利用載波相位來進行定位；其次，載波相位觀測量的量測可以非常精準，而電碼相位觀測量的量測則較為粗糙。舉例而言，市面上一些較高階的GPS接收機，若不考慮其他誤差的影響，載波相位量測誤差的標準差大約可以精確至$\sigma ({\varepsilon _\phi }) \cong 0.025{\rm{ cycle}}$(5毫米)，而電碼相位觀測量誤差的標準差大約為$\sigma ({\varepsilon _\phi }) \cong$0.5 m。

GPS載波相位觀測量大多是採用相對定位 (relative positioning) 的方式，利用在兩個不同地點的GPS相位觀測量，計算兩地之間的相對位置向量 (relative position vector)，在測量學上稱此相對位置向量為基線向量 (baseline vector)，或者簡稱為基線。利用GPS進行精密的相對定位，需要在兩個不同點上同時接收載波相位觀測量，天線在參考位置 (座標為已知) 上的接收機稱為參考接收機 (reference receiver) 或者參考站，而第二架接收機天線位置係在待測的點上，稱為移動接收機 (mobile receiver or rover) 或者移動站。利用載波相位相對定位需要一個所謂”初始化” (initialization) 的過程，亦即求解整數週波未定值，若是在初始化以及定位的過程中，移動站以及參考站皆為靜止狀態，稱為靜態測量 (static survey)。事實上，在相對定位過程中，一旦求得整數未定值，移動站即可自由移動，因此靜態測量實際上僅要求在初始化的過程中移動站是靜止的。而在初始化的過程中，若移動站可以自由移動則稱為動態測量 (kinematic survey)。早期 (1980年代)，動態初始化的方法尚未發現，靜態GPS測量通常需要接收機處於靜止的狀態，接收大約一小時左右的載波相位觀測量之後，等待衛星幾何分佈有明顯的改變，以求解整數未定值。在1993年，學者Talbot最早提出動態初始化的方法[5]，可以動態地進行相對定位，使得大地測量更有效率，此種方法稱為即時動態定位 (real-time kinematic positioning，簡寫為RTK)，在RTK模式下，參考站透過適當的數據鏈，將載波相位觀測量即時地傳送至移動站。RTK技術的關鍵是移動站必須有能力可以在運動的狀況下，估測整數未定值，此種方式稱為On-The-Fly (OTF) Initialization。

(待續)

References:

1. C.C. Counselman III, I.I. Shapiro, R.L. Greenspan, and D.B. Cox, Jr., “Backpack VLBI Terminal with Subcentimeter Capability,” Proc. Radio Interferometric Techniques for Geodesy, NASA Conference Publication, Vol. 2115, pp. 409-413, 1979.
2. C.C. Counselman III and Sergei Gourevitch, “Miniature Interferometer Terminals for Earth Surveying: Ambiguity and Multipath with Global Positioning System,” IEEE Transactions on Geosciences and Remote Sensing, Vol. GE-19, No. 4, pp. 244-252, 1981.
3. Patrick Y.C. Hwang, “Kinematic GPS for Differential Positioning: Resolving Integer Ambiguities on the Fly,” Navigaiton, Vol. 28, No. 1, pp. 1-15, 1991.
4. Benjamin W. Remondi, “Centimeter-Level Surveys in Seconds with GPS Carrier Phase: Initial Results,” Navigaiton, Vol. 32, No. 4, pp. 386-400.
5. Nicholas Talbot, “Centimeters in the Field, A User’s Perspective of Real-Time Kinematic Positioning in a Production Environment,” Proc. ION GPS-93, pp. 1049-1057, 1993.

第3章 GPS座標與時間系統 (II)

第3.2節 全球座標系統
3.2 全球座標系統
座標系統的定義，在導航系統中佔有重要的角色，座標系統就其涵蓋範圍來區分，應可大致分為本地座標系統 (Local-Level Coordinate System) 以及全球座標系統 (Global Coordinate Systems)，在日常生活中，一般大眾並不太需要用到所謂的全球座標系統 (例如我與朋友相約，我不會跟他說我們約在東經121度52分，北緯24度31分的地方見面)。但是，作為一個全球通用的導航系統，GPS系統必須定義一個全世界任何地方均可適用的座標系統，也就是全球座標系統。全球座標系統必須可以用來描述地球上任意一點的三度空間座標，在GPS導航系統中，常用的全球座標系有兩種：地心固連 (或者稱為地心地固Earth-Centered, Earth-Fixed; ECEF) 座標系統以及地球慣性 (Earth-Centered, Inertial; ECI) 座標系統。全球座標系統是一種卡氏座標系統，由我們在上一節中的討論可知，定義一個卡氏座標系統最主要的元素包括原點以及三個互相垂直的座標軸。所以，我們在定義地心固連座標系統以及地球慣性座標系統時，實際上就是在指定其原點所在位置，以及三個座標軸的指向。

3.2.1 地球座標系統與慣性座標系統
在定義地球座標系統時，經常會用到所謂參考子午線以及參考橢球的概念，因此我們先將其定義整理如下：

定義3.1 (參考橢球) 給定一平面上之橢圓，其半長軸a為地球平均半長軸而且與x軸平行，半短軸b為地球平均半短軸而且與z軸平行，若將此橢圓延地球自轉軸旋轉180o，所形成之球體即稱為大地參考橢球 (Reference Ellipsoid) 或者簡稱為參考橢球 (見圖3.6)。參考橢球之球心位於地球的地理中心，其與xz平面相切之切面稱為參考子午線(Reference Meridian)，與xy平面相切之切面稱之為平均赤道面 (Mean Equator) 或者簡稱為赤道面。

在上面的定義中，除了z軸定義為地球的自轉軸之外，x軸與y軸尚未定義，事實上，根據x軸定義的不同，即產生不同的大地座標系統，常用的定義方式有兩種。首先我們考慮第一種定義方式，假設地球的質量中心為座標原點，地球的自轉軸為z軸，而x軸通過赤道面與格林威治平均子午線 (Mean Greenwich Meridian) 的交點，y軸則由x軸與z軸以右手定則決定，則此座標系即稱為地心固連坐標系統。這個定義方式最大的問題是地球的自轉軸實際上並不是固定的，而是進行一種近乎週期性的運動，稱之為極軸運動效應或者極動 (polar motion)，因為赤道面與自轉軸垂直，因此赤道面也會跟著移動，而一旦赤道面移動，經緯線的定義也會跟著變動。圖3.7所示為2000年至2005年間之瞬時極軸位置圖。

為解決這個問題，可以將z軸重新定義為地球的平均自轉軸，大地測量學者利用西元1900年至1905年間之極軸運動效應，定義所謂的協議地極 (Conventional Terrestrial Pole; CTP)，以此點定義地球的平均自轉軸。而與此自轉軸垂直之赤道面稱為CTP赤道面。

定義3.2 (協議地球座標系統) 考慮如圖3.8所示之參考橢球與卡氏直角座標系統，其原點位於地球的質量中心。若其z軸通過CTP，x軸通過CTP赤道面與參考子午線的交點，y軸由x軸與z軸以右手定則決定構成一直角座標系統，則此座標系統稱之為協議地球座標系統 (Conventional Terrestrial Reference System; CTRS)，是一種地心地固座標系統。此時，其x, y, z軸分別以$x_T, y_T, z_T$表示。

上述的定義屬於”抽象”的說法，舉例而言，實際上我們並無法”看”到地球的質心所在位置，而且參考子午線、赤道面等，也都是虛擬的線以及平面，必須透過大地量測的方法，才能將一個大地座標系統”實現” (realization) 出來。例如WGS84即為CTRS的一種實現，我們將在後面的章節介紹WGS84參考座標系統。在這裡必須先提出一項說明，由上述的大地量測方法所實現出來的座標系統，其零度經度線不再是通過格林威治天文觀測站那條著名的銅線，而是利用位於全球各地的天文觀測站，經過協調之後統計上的零度線，亦即所謂的格林威治平均子午線。

CTRS座標系統的x軸指向格林威治子午線，因此CTRS座標系統會隨著地球自轉而跟著轉動，所以CTRS並非慣性座標系統。此種座標系統適合用於描述地表附近物體的位置，但並不適合用來描述與分析衛星軌道運動。衛星軌道運動為一種向心力運動 (Central Force Motion)，其運動軌跡可以由牛頓運動定律來決定，而描述牛頓運動定律使用慣性座標系統是最佳的選擇，因為在此座標系統之下，牛頓運動定律有最簡潔的形式。接下來我們就來介紹另一種全球座標系統：地球慣性 (Earth-Centered, Inertial; ECI)座標系統。

定義3.3 地球繞太陽公轉之軌道面稱為黃道面 (Ecliptic)，而黃道面與赤道面交線所指方向即為春分點 (Vernal Equinox, 見圖3.9)。

定義3.4 (協議慣性座標系統) 考慮如圖3.10 協議慣性座標系統 (CIRS)所示之參考橢球與卡氏直角座標系統，其原點位於地球的質量中心。若其z軸通過地球的瞬時自轉軸(Celestial Ephemeris Pole; CEP)，x軸通過春分點，y軸由x軸與z軸以右手定則決定，則此座標系統稱之為協議慣性座標系統 (Conventional Inertial Reference System; CIRS)。此時，其x, y, z軸分別以$x_I, y_I, z_I$表示。

嚴格說來，上述定義並非一個真正的慣性座標系統，因為地球的質心會繞著太陽進行向心運動，即所謂的公轉，不過在實用上，在短時間之內CIRS可以視為慣性座標系統。慣性座標系統也可以稱為太空固定 (space-fixed) 系統，因為其x軸指向在太空中固定不動的恆星。除了質心會移動之外，地球的自轉軸相對於遠處的恆星而言也並非固定不動的，地球自轉軸的運動主要是由兩種不同的週期運動所構成，即所謂的進動 (precession) 與章動 (nutation)，這兩種運動是由於太陽與月亮對地球的引力所造成的。進動就是地球的自轉軸本身，會進行繞圓圈的運動 (其週期大約26,000年)，就像陀螺一樣，章動則是自轉軸會像點頭一樣，前後擺動 (週期大約18.6年)。如果地球是正圓形而且為均勻物質 (密度固定)，就不會產生進動與章動的現象，不過幸好科學家對於這兩種運動已經有充分的瞭解，可以準確的計算任何時間點 (epoch) 地球的進動與章動。

為了描述地球在太空中的方位 (姿態) 以及CIRS與CTRS兩種座標系統之間的相對關係，除了進動與章動之外，還需考慮另外兩種因素，首先是前面提到的地球的極軸運動效應，另外一個因素則是地球的自轉速率，因為地球的自轉速率實際上並不是固定不變的。綜合上述的幾項因素，我們需要五個參數來描述CTRS與CIRS之間的關係，以及地球的方位，這五個參數稱為EOP (Earth Orientation Parameter)，整理如下：

• 兩個角度($x_p$以及$y_p$) 用以描述極軸運動效應座標，亦即地球的瞬時自轉軸與CTP之間的相對位置，如圖3.11所示，其中$(x_p,y_p)$代表極軸運動效應座標 (Polar Motion Coordinate, 以CTRS座標表示)。
• 利用兩個角度來表示地球自轉軸在太空中的指向，而且必須考慮進動與章動對於地球自轉軸的影響。
• 以$\theta$角代表格林威治平均子午線與春分點之夾角，用來描述地球的自轉運動。$\theta$角的定義必須同時考慮地球非均勻自轉率的影響 (見圖3.11)。

圖3.11 CTRS與CIRS座標關係圖

地球的極軸運動效應與進動及章動不一樣的地方是到目前為止極軸運動效應仍無法精確預測，因此，CEP相對於CTP之位移  僅能利用實驗的數據來判斷當時的極軸運動座標值。$\theta$角代表格林威治平均子午線與春分點指向之夾角，也稱為格林威治視恆星時 (Greenwich Apparent Sidereal Time; GAST)。一旦地球的進動與章動可以精確的計算，我們僅需要極軸運動座標以及地球自轉角度三個參數  來描述CIRS以及CTRS兩個座標系統的相對關係。

由本章節附錄 (見第3.6節) 的討論中可知，兩個有共同原點的不同卡氏直角座標系統，可以經由一系列的旋轉矩陣 (rotational matrix，或者稱為方向餘弦矩陣Directional Cosine Matrix; DCM) 運算之後，使得兩個座標系統重疊，換句話說不同卡氏直角座標系統之間的相對關係可以由方向餘弦矩陣來表示。我們可以利用圖3.11來說明如何將慣性座標$(x_I,y_I,z_I)$轉換至ECEF座標$(x_T,y_T,z_T)$，首先將$z_I$軸固定，順時針旋轉$\theta$角，接著相對於新的$x_I$軸旋轉$-y_p$角，最後再固定新的$y_I$軸旋轉$-x_p$角，由此三次旋轉可得知CIRS系統與CTRS系統之間的座標轉換矩陣可以表示如下：

\[{\bf{x}}_T = {\bf{R}}_{2I}( - {x_p}){\bf{R}}_{1I}( - {y_p}){\bf{R}}_{3I}(\theta ){\bf{x}}_I.\]

其中$x_T$與$x_I$分別代表以CTRS與CIRS座標系統表示之位置向量，${\bf R}_{1I}$、${\bf R}_{2I}$以及${\bf R}_{3I}$分別為相對於$x_I$、$y_I$以及$z_I$軸旋轉之方向餘弦矩陣。在本書中，粗體的大寫英文字母代表矩陣，粗體小寫英文字母則代表向量。

3.2.2 大地座標系統、大地水準面以及大地基準
數學上，卡氏直角座標系統有很多好處，例如牛頓運動定律以及馬克士威方程式 (Maxwell’s Equations) 在此座標系統之下，公式較為簡潔。但是在日常生活上，卡氏直角座標系統使用起來卻嫌累贅。想像GPS接收機告訴你現在所在位置為 (1,510,900, -4,463,457, 4,293,001) (以公尺為單位)，大部分的人大概知道是在北半球 (因為z > 0)，如果你記得三角函數表，或者你的算術很好，可以用心算得知你現在位於中緯度 ($40^\circ-45^\circ$) 地區。不過，如果要判斷所在地點的高度，可能就沒有那麼容易了！在地球表面上使用直角座標的另一個缺點是就算你的高度僅有些微改變 (例如爬了一層樓)，三個座標值都會改變。因為地表接近球型，在日常生活中採用曲線座標系統 (經度longitude、緯度latitude、高度height) 反而更為便利，也就是所謂的大地座標系統 (Geodetic Coordinate System)。

大地座標系統的概念是先以一個橢球模型 (見定義3.1) 做為地球表面的近似模型，橢球表面上的任意一點僅需使用兩個角度來描述 – 經度與緯度，高度則是以相對於橢球表面的高度 (即所謂的海拔altitude) 來表示。地球的表面並不規則，甚至會隨著時間改變，因此要找到一個適當的橢球模型並不容易。早期的模型是一個球型，如圖3.12所示，由此定義出來的經度與緯度稱為地心經度 (Geocentric Longitude) 與地心緯度 (Geocentric Latitude)。一直要到牛頓發現萬有引力定律之後，人類才知道地球實際上是一個扁的橢球形 (Oblate Ellipsoid)，如圖3.13所示 (橢球的定義請參考定義3.1)，根據測量的結果，地球的長軸與短軸大約相差20公里。

由上面的敘述可知，大地座標系統除了一組ECEF座標軸以及原點之外，還需要一個參考橢球模型，而為了簡化模型複雜度，可令橢球球心與ECEF原點重疊，而橢球的旋轉軸與ECEF座標的z軸重合。有了橢球的旋轉軸與球心之後，由定義3.1可知，我們還需要兩個參數才能完整描述一個參考橢球，亦即半長軸($a$)以及半短軸($b$)，通常也以半長軸以及離心率($e$)來表示。離心率與半長軸以及半短軸的關係如下：
\[{e^2} = \frac{{({a^2} - {b^2})}}{{{a^2}}}.\]
大地測量學家則比較喜歡用扁率 (flattening, $f$) 來描述橢球，扁率的定義如下：
\[f = \frac{{(a - b)}}{a}.\]
由上面兩個式子可以得到離心率與扁率的關係，如下列方程式所示：
\[{e^2} = 2f - {f^2}.\]
接著，我們可以來定義大地座標(也稱為地理座標Geographic Coordinate，或者橢球座標Ellipsoidal Coordinate)。

定義3.5 (大地座標系統) 考慮圖3.14所示之直角座標與參考橢球，空間中的一點P其卡氏直角座標為 (x, y, z)，則其大地座標定義如下：

• 大地緯度 (geodetic latitude, φ)：P點所在的子午面上，赤道面 (x-y平面) 與P點在橢球表面的垂線之間的夾角，赤道北方緯度為正，南方為負。
• 大地經度 (geodetic longitude, λ)：在赤道面上，參考子午面以及P點所在子午面之間的夾角 (兩面角)，參考子午面 (零度子午面) 的東方經度為正。
• 大地高程 (geodetic height, h)：P點到橢球表面之間的垂直距離。

在圖3.14中，φ’ 代表地心緯度，係指$\overline {{\bf{OP}}}$與赤道面之夾角 (O為原點)。因為地球為橢球形，所以。若地表為正球形，則。但是因為實際上地球為非常接近球形的橢圓形，故一般說來$\phi  \approx \phi '$。在日常生活中，我們所謂的經緯度係指大地經緯度。將橢球表面上，所有緯度相同 (φ = constant) 的點連接其來形成的一條封閉曲線，稱為平行圈 (parallel)，所有平行圈都是圓形的。同理，將橢球上面所有經度相同 (λ = constant) 的點連接起來形成的封閉曲線，稱為子午圈 (meridian)，所有子午圈都是橢圓形。大地高程也稱為橢球高 (ellipsoidal height)，也就是P點相對於橢球面的高度，因為大地參考橢球面為一個假想的曲面，所以橢球高並沒有實際上的物理意義。因此大地測量學家定義另一種曲面，稱為大地水準面 (Geoid)。

定義3.6 考慮圖3.15。鉛垂線 (Plumb Line)、等位面 (Equipotential Surface) 以及大地水準面 (Geoid) 可分別定義如下：

• 等位面：連接所有重力位能相同的點所形成之曲面。
• 鉛垂線：與等位面垂直的線稱為鉛垂線，同時也是地心引力的方向。
• 大地水準面：與全球平均水平面 (global mean sea level) 最接近之等位面。

說明：

1. 在上述定義中，所謂”最接近”係指以最小平方法 (least squares method) 進行曲線契合 (curve fitting) 所得到之最佳解 (optimal solution)。
2. 鉛垂線也稱為垂直線 (Vertical)，與橢球面垂直的線則稱為法線 (Normal)。

定義3.7 參考圖3.16，在地形上的任意一點P其橢球高 (Ellipsoidal Height)、正高 (Orthometric Height) 與大地水準面高 (Geoidal Height) 之定義分別為：

• 橢球高：P點到參考橢球面之法線長度，亦即圖中的h。
• 大地水準面高：橢球面以及大地水準面之法線距離，亦即圖中的N。
• 正高：地表與大地水準面之鉛垂線長，以大寫H表示。

圖3‑17 大地水準面、大地水準面高以及垂線偏角 (Deflection of the Vertical)

說明：

1. 對於地表的任意一點P，通過P點的法線與鉛垂線的夾角稱為垂線偏角 (Deflection of the Vertical，見圖3.17)。
2. 一般說來，橢球高、正高以及大地水準面高有下列關係：\[h=H+N. \]上式實際上是一個近似關係 (approximate relation)，不過在大部分情況下其精確度已經足夠。
3. 橢球高以及大地水準面高均以橢球面為參考面，橢球面上方為正，下方為負，亦即在圖3.16中，N值為負。

在做地圖投影 (Map Projection) 時，比較方便的方式就是將3D的地形投影到參考橢球面上，地圖投影的方式已經超過本書涵蓋範圍，在此不詳細介紹。在這裡我們特別要注意的是，如果是區域性的地圖投影，使用全球參考橢球為基準，實際上並不合適，使用區域基準面 (Regional Datum) 可能更為精確。大多數的國家均會發展自己的區域座標系統以及參考橢球，台灣目前採用的大地基準，稱為TWD 97系統，將在後面章節介紹。大地基準與區域基準的差異如圖3.18所示。台灣的大地水準面高輪廓圖如圖3.19所示。

參考文獻

1. P. Misra and P. Enge, Global Positioning System: Signals, Measurements, and Performance, Ganga-Jamuna Press, 2001.
2. 王和盛, GPS接收機原理概論, unpublished manuscript.

Nonbinary m-Sequences

m-sequences 在衛星導航系統中扮演重要的腳色，測距用的虛擬雜訊碼 (Pseudo-Random Code)基本上就是一種m-sequence。有關GPS訊號中PRN碼的設計可參考下面這個鏈結中的文件。

GPS訊號結構

一般的m-sequence大多為二進制(binary)，但也有非二進制的m-sequence。本篇主要介紹非二進制的m-sequence。欲生成一個非二進制的m-sequence，需要先建構一個extension field $GF(p^m)$，其中$p$是一個質數。在工程上最常用的extension field是$GF(2^m)$，可參考上述文章，有詳細介紹如何生成$GF(2^m)$。我們以一個實際的例子來說明，考慮一個二階的primitive polynomial $X^2+X+1$，則以此多項式生成之$GF(2^2)$ extension field包含四個元素，以冪次表示為$\{0, 1, \alpha, \alpha^2\}$，以多項式表示則為$\{0, 1, \alpha, \alpha^2 = 1 + \alpha\}$。也可以利用向量來表示，此時$GF(2^2) = \{(0,0), (0,1), (1,0), (1,1)\}$。接下來，我們考慮在$GF(2^2)$上的多項式，同樣地，為了方便說明，我們考慮一個二階的primitive polynomial $g(X)=g_0+g_1X+g_2X^2 = X^2+X+\alpha$。注意，因為我們現在考慮的是$GF(2^2)$上的多項式，因此$g(X)$的係數可以是$0, 1, \alpha,$ 或者$\alpha^2$。在本例中$g_0=\alpha, g_1 = 1, g_2=1$。接下來我們可以利用$g(X)$來產生一組非二進制m-sequence。若令$g(X)=0$，移項之後可得$X^2= X + \alpha$。若將$X$視為電路中的移位暫存器，$X^2$表示經過兩個移位暫存器(有關以移位暫存器來產生一個序列，同樣可參考上述連結中的文章)，則由上式可得到一個遞迴關係：$a_r = a_{r-1} + \alpha a_{r-2}$。若令$a_0=a_1=1$，由此遞迴關係可得到一個長度為$2^4-1 = 15$的非二進制m-sequence如下：
\[a_0a_1a_2a_3... = 11\alpha^210\alpha\alpha 1\alpha 0\alpha^2\alpha^2\alpha\alpha^2 0 \mbox{(repeat)} \]
但是，若想要利用數位電路來實現上述m-sequence，要如何產生$\alpha$這個符號，此時就需要借用$GF(2^2)$的向量表示，也就是說$\{0, 1, \alpha, \alpha^2\}$這四個符號都需要兩個位元來表示：$"0"\rightarrow 00$, $"1"\rightarrow 01$, $"\alpha"\rightarrow 10$, $"\alpha^2"\rightarrow 11$。如此一來，上述的非二進制m-sequence可以寫為：
\[01,01,11,01,00,10,10,01,10,00,11,11,10,11,00\]
解碼時必須兩個位元為一組進行解碼。

References
1. 王和盛, GPS接收機原理概論, unpublished manuscript.
2. S. Lin and D. J. Constello, Jr., Error Control Coding: Fundamentals and Applications, 2nd Edition, Prentice-Hall.

隨機程序 Random Processes

在介紹隨機程序(或稱隨機過程)之前，先簡單地回顧一下機率理論。若一個實驗的結果(outcome)無法事先得知，即可稱為隨機實驗(Random Experiment)，這是在機率理論的一個重要基本概念。在一個隨機實驗中，所有可能的結果形成的集合(set)，稱為該實驗的樣本空間，並以$\Omega$表示。樣本空間的任一個子集合(subset)，稱為一個事件(event)，記為$E$。若實驗的結果為某一事件$E$的一個元素，則稱發生該事件。對樣本空間中的每一個事件$E$，我們定義一個集合函數$P(E)$，且該函數滿足下列三個公設(axioms)：

• Axiom 1: $0\le P(E)\le 1$.
• Axiom 2: $P(\Omega) = 1$.
• Axiom 3: 對任意一組事件序列$\{E_i\} = E_1, E_2, E_3,...$ (可以有無限多個)，若對任意的$i\ne j$，$E_i\cap E_j = \emptyset$，亦即$\{E_i\}$為互斥(mutually exclusive)，則\[ P\left(\bigcup_{i=1}^\infty E_i\right) = \sum_{i=1}^\infty P(E_i).\]

若這樣的函數存在，則我們稱$P(E)$為事件$E$發生的機率。令$\mathcal{F}$表示所有可能事件所形成的集合($\mathcal{F}$嚴格說來，必須是一個$\sigma$-field，或者稱為$\sigma$-algebra，不過在這邊暫時可以不管這些，那為何要提這個咧？因為這樣看起來好像比較專業。)，則$(\Omega, \mathcal{F}, P)$這三個物件構成一個機率空間(Probability Space)。這種機率的定義方式，稱為公設式的定義(axiomatic approach)，是數學上最常用的定義方式。公設有人喜歡有人不喜歡，像是在台北市買房子，有的還規定公設至少30%，真是他媽的坑爹啊！數學家則是非常喜歡公設，機率的公設是由俄羅斯的數學家Kolmogorov提出來的，他告訴我們，滿足這三個公設就是機率，也就是說這個世界上只有機率可以同時滿足這三個性質，不需要計算事件發生的頻率，也不需要假設你的實驗公不公平。這種公設式的定義，在數學上隨處可見，有時候不得不說：這些數學家真是他媽的有病啊！不然怎麼在一堆雜亂無章的數學式中，找出這兩、三個規則，並且將它定為不可違背的公設。下面列出機率函數的幾個簡單的性質，利用三個基本公設即可證明。

1. If $E\subset F$, then $P(E)\le P(F)$.
2. $P(E^c) = 1 - P(E)$, where $E^c$ is the complement of $E$.
3. $P(\bigcup_{i=1}^n E_i) = \sum^n_{i=1} P(E_i)$ when the $E_i$ are mutually exclusive.
4. $P(\bigcup_{i=1}^n E_i) \le \sum^n_{i=1} P(E_i)$.

接下來介紹隨機變數(Random Variables)。考慮一個樣本空間$\Omega$，令$\omega$表示一次隨機實驗的結果，隨機變數$X: \Omega\rightarrow\mathbb{R}$是一個函數，指定一個唯一的數值$\xi$給該次實驗結果，亦即$X(\omega) = \xi$(參考圖1)。若$X$的值域為$\mathbb{R}$則稱為連續隨機變數，若$X$的值域為$\mathbb{Z}$或$\mathbb{N}$，則稱為離散隨機變數。事實上，隨機變數還需要滿足另一個條件，如下：\[ P\{X\in A\}=P(X^{-1}(A)). \] 其中$A$為$\mathbb{R}$的一個子集合，$X^{-1}(A)$表示$A$的逆映像(inverse image)，簡單來講就是對任意的$A\subset \mathbb{R}$，其inverse image必須是一個事件，且發生$A$的機率與發生$X^{-1}(A)$的機率相同。同樣地，這個性質可以先暫時不管他，我只是想告訴大家，我很專業！(個屁)

圖1：隨機變數

隨機變數的概念，可以說是近代機率理論的重要基石，有些讀者和同學(包括我也是)一開始可能會覺得隨機變數的概念不容易接受，而且似乎有點多餘，不過是要將每個事件，指定一個數值，為何要特別定義一個函數咧？初學時，其實可以只要先記得隨機變數是一個函數這樣的概念就好了，它的很多好處需要慢慢體會，我自己的認知，隨機變數的概念可以簡化許多機率的問題，任何無法預知(random or stochastic)的東西，經過隨機變數映射之後(也就是隨機試驗的一個"實現", realization)，都變成是"確定性"(deterministic)的東西，一但是確定的東西，我們就可以對他進行微分、積分以及其他複雜的數學運算。這個概念在討論隨機程序時，更能顯示出他的好處。

令$\Omega$為樣本空間(Sample Space)，$T$為實數的一個子集合，代表一個時間區段(time interval)。則隨機程序(random process or stochastic process)可以定義為從product space $T\times\Omega$(亦即兩個變數的函數，其中第一個變數為$T$中的元素，第二個變數為$\Omega$中的元素)映射至$\mathbb{R}^n$的一個函數，使得對任意的$t_j\in T$，$\mathbf{X}(t_j,\cdot)$是一個隨機變數($n=1$時是純量隨機變數，$n>1$為多變量隨機變數)。反之，若固定第二個變數，則對任意的$\omega_i\in\Omega$，$\mathbf{X}(\cdot,\omega_i):=\mathbf{X}(\cdot)$是一個時間函數，稱為樣本函數(sample function)。對固定的 $t_j\in T$，$\omega_i\in\Omega$，則得到$\mathbf{X} (t_j,\omega_i)$是一個確定的向量。

圖2: 隨機程序

Gauss-Markov Process

考慮一個一階的線性常微分方程式：
\[\dot{x}(t)+ax(t)=0,\;\; x(0)\;\; \mbox{is given.}\]
其解為
\[x(t)=e^{-at}x(0). \]
假設採樣時間(sampling time)為$\Delta T$，且令$t_k:=k\Delta T$，則不難證明，由上述方程式可得
\[ x(t_{k+1}) = e^{-a\Delta T}x(t_k). \]
上式可視為一個離散時間的線性系統，其時間常數(time-constant)為$1/a$。若令此系統的輸入$w(t_k)$，為一個zero-mean Gaussian random sequence (or process)，其變異量為$\sigma^2$，此時上式可改寫為：
\[x(t_{k+1}) = e^{-a\Delta T}x(t_k)+w(t_k). \]
則Gauss-Markov定理告訴我們，$\{x(t_k)\}_{k=0}^{k=\infty}$必為一個Gauss-Markov process，且$x(t_k)$的autocorrelation function為
\[R(\tau)=\sigma^2e^{-a\vert\tau\vert}. \]

可以利用Gauss-Markov Process來描述GPS衛星誤差模型。上圖中Error Parameter就是指$x(t)$，Standard Deviation = $\sigma$，Time = 時間常數。舉例而言，第二列的意思就是說將Ephemeris以$\sigma=3, a=1/1800$的Gauss-Markov Process來表示，其autocorrelation function $R(\tau) = 9e^{-\vert\tau\vert / 1800}$。

Gram-Schmidt正交化與傅立葉轉換之間的關係

正交化是內積空間上的核心觀念，假設有一組獨立的向量集合，則可利用Gram-Schmidt正交化過程，從這組獨立集，造出一組正交集(orthogonal set)以及單範正交集(orthonormal set)。由此方法可推出矩陣的兩種類型QR分解，Gram-Schmidt正交化概念，與傅立葉轉換兩者之間有密不可分的關係。

定理1. (Fourier Coefficient)
設$S=\{\mathbf{v}_1,\mathbf{v}_2,..., \mathbf{v}_k\}$為內積空間$\mathcal{V}$之非零向量正交集，且$\mathbf{v}\in\mbox{span}(S)$。若$\mathbf{v}=\sum^k_{i=1}\alpha_i\mathbf{v}_i$，則
\[\alpha_j = \frac{<\mathbf{v},\mathbf{v}_j>}{\Vert\mathbf{v}_j\Vert^2},\;\;\;\;\mbox{for all }j=1,2,..., k.\]
$\alpha_i$稱為傅立葉係數(Fourier Coefficient)。

定理2.(Gram-Schmidt Orthogonalization Process)
已知$S=\{\mathbf{v}_1,\mathbf{v}_2,..., \mathbf{v}_k\}$為內積空間$\mathcal{V}$之獨立集，則令
\[\begin{array}{rcl}
\mathbf{u}_1 & = & \mathbf{v}_1, \\
\mathbf{u}_2 & = & \mathbf{v}_2 - \frac{<\mathbf{v}_2,\mathbf{u}_1>}{<\mathbf{u}_1,\mathbf{u}_1>}\mathbf{u}_1, \\
\vdots & & \vdots \\
\mathbf{u}_k & = & \mathbf{v}_k - \sum^{k-1}_{i=1} \frac{<\mathbf{v}_k,\mathbf{u}_i>}{<\mathbf{u}_i,\mathbf{u}_i>}\mathbf{u}_i,
\end{array}\]
則$S'=\{\mathbf{u}_1,\mathbf{u}_2,..., \mathbf{u}_k\}$形成非零向量正交集。

假設$S=\{\mathbf{v}_1,\mathbf{v}_2,..., \mathbf{v}_n\}$為內積空間$\mathcal{V}$的一組基底，則可以利用定理2，造出一組正交基底$S'=\{\mathbf{u}_1,\mathbf{u}_2,..., \mathbf{u}_n\}$。依據基底的定義，對於$\mathcal{V}$上的任意向量$\mathbf{v}$，皆可以表示成$\mathbf{v}=\sum^n_{i=1}\alpha_i\mathbf{u}_i$，且
\[\alpha_j = \frac{<\mathbf{v},\mathbf{u}_j>}{\Vert\mathbf{u}_j\Vert^2},\;\;\;\;\mbox{for all }j=1,2,..., n.\]
$\alpha_j$稱為傅立葉係數，實際上就是向量$\mathbf{v}$在基底$\mathbf{u}_j$上的分量大小。現在，我們將這個概念，推廣到訊號空間上。在繼續之前，需要先說明一下向量空間維度的概念。$\mathbb{R}^n$空間的維度為$n$，所代表的意思就是說我們只需要$n$個線性獨立的向量即可完整描述$\mathbb{R}^n$空間，也就是說$\mathbb{R}^n$空間中的任一個向量，都可以用這$n$個向量的線性組合來表示。那麼，訊號空間的維度是多少呢？從線性代數的理論可以證明，訊號空間的維度是無窮大，因此我們需要無限多個基底才能完整描述訊號空間中的任意"向量"(意即信號)。若基底的個數為有限，則僅能以近似的方式表示訊號空間中的向量。

假設$\phi_1(t), \phi_2(t),..., \phi_N(t)$為訊號空間中的一組正交向量，

Impulse Signal Generation - MATLAB

Impulse signal generation—MATLAB

When defining the impulse or $\delta (t)$ signal, the shape of the signal used to do so is not important. Whether we use the rectangular pulse we considered in this chapter or another pulse, or even a signal that is not a pulse, in the limit we obtain the same impulse signal. Consider the following cases:
(a) The triangular pulse,
\[\Lambda_{\Delta}(t) = \frac{1}{\Delta}\left(1-\left\vert \frac{t}{\Delta}\right\vert\right)\left(u(t+\Delta) - u(t-\Delta)\right). \]
Carefully plot it, compute its area, and find its limit as $\Delta\rightarrow 0$. What do you obtain in the limit? Explain.
(b) Consider the signal
\[S_\Delta (t) = \frac{\sin (\pi t/\Delta)}{\pi t}. \]
Use the properties of the sinc signal $S(t) = \sin (\pi t)/(\pi t)$ to express $S_\Delta (t)$ in terms of $S(t)$. Then find its area, and the limit as $\Delta\rightarrow 0$. Use symbolic MATLAB to show that for decreasing values of $\Delta$ the $S_\Delta (t)$ becomes like the impulse signal.

Solution: MATLAB Script
(a)
% Pr. 1.7
clear all; clf
% part (a)
delta=0.1;
t=[-delta:0.05:delta];N=length(t);
lambda=zeros(1,N);
figure(5)
for k=1:6,
lambda=(1-abs(t/delta))/delta;
delta=delta/2;
plot(t,lambda);xlabel(’t’)
axis([-0.1 0.1 0 330]);grid
hold on
pause(0.5)
end
grid
hold off

(b)
% part (b)
syms S t
delta=1;
figure(6)
for k=1:4,
delta=delta/k;
S=(1/delta)*sinc(t/delta);
ezplot(S,[-2 2])
axis([-2 2 -8 30])
hold on
I=subs(int(S,t,-100*delta, 100*delta)) % area under sinc
pause(0.5)
end
grid;xlabel(’t’)
hold off

向量分析 (III)

1.3 向量積分

1.3.1 線、面以及體積分

In electrodynamics we encounter several different kinds of integrals, among which the most important are line (or path) integrals, surface integrals (or flux), and volume integrals.
(a) Line Integrals. A line integral is an expression of the form
\[\int^\mathbf{b}_{\mathbf{a}\mathcal{P}}\mathbf{v}\cdot \mathrm{d}\mathbf{l}, \]
where $\mathbf{v}$ is a vector function, $d\mathbf{l}$ is the infinitesimal displacement vector, and the integral is to be carried out along a prescribed path $\mathcal{P}$ from point $\mathbf{a}$ to point $\mathbf{b}$. If the path in question forms a closed loop (that is, if $\mathbf{b}=\mathbf{a}$), we shall put a circle on the integral sign:
\[\oint\mathbf{v}\cdot d\mathbf{l}. \]
To a physicist, the most familiar example of a line integral is the work done by a force $\mathbf{F}$:
\[\int\mathbf{F}\cdot d\mathbf{l}. \]
Ordinarily, the value of a line integral depends critically on the particular path taken from $\mathbf{a}$ to $\mathbf{b}$, but there is an important special class of vector functions for which the line integral is independent of the path, and is
determined entirely by the end points. A force that has this property is called conservative.

$$$$
∑i=1Nai

訊號與系統 (II)

1.2 幾種有用的信號操作方式
1.2.1 時間位移
考慮一信號$g(t)$，以及該信號延遲$T$秒的版本，以$\phi (t)$表示。則
\[\phi (t+T) = g(t)\;\;\;\;\;\mbox{或者} \;\;\;\;\; \phi(t) = g(t-T). \]
因此(見圖1.2)，若要將一個信號在時間上平移$T$，我們將$t$以$t+T$取代，亦即$g(t-T)$表示將信號$g(t)$平移$T$秒。若$T$是正的，代表係往右平移(延後)，若$T$是負的，表示信號係往左平移(超前)。

圖1.2 信號的時間平移

1.2.2 時間刻度變換

第3章 GPS座標與時間系統 (I)

GPS衛星繞著地球進行夾角約55度的軌道運動，其位置會隨著時間而改變。從幾何的觀點來看，軌道座標系統最適合用來分析衛星的運動狀態。另一方面，GPS定位或是量測大多數是在地表進行，因此使用地心座標系統較為方便，而且因為地表接近圓形，最適合用來描述地表上任一點位置的方式是使用極座標系統，而非直角座標系統。GPS定位需要量測衛星所在位置 (軌道座標系統、隨時間改變) 以及接收機所在位置 (極座標系統) 之間的距離，以計算出使用者所在位置，在此過程中，需要處理不同座標系統以及不同的參考時間資料，因此如何定義與選擇一個適當的座標軸以及時間系統，是設計導航系統的重要課題。

座標與時間系統與天體力學 (celestial mechanics) 之間的密切關係，可以回溯至數百年前，在十七與十八世紀時，許多偉大的科學家，例如：伽利略 (Galileo)、克卜勒 (Kepler)、牛頓 (Newton)、高斯 (Gauss)、惠更斯 (Huygens) 以及尤拉 (Euler) 等人，各有許多偉大的貢獻，其中尤其是牛頓在天體軌道動力學的貢獻更是眾所公認的。古典力學是描述具有質量的物體受力時，其運動方式的變化以及相互間的規則。在牛頓 (Isaac Newton, 1642-1727) 之前，描述物體運動的基本觀念與數學方法皆已大致建立，其中的代表性人物是克卜勒 (Johannes Kepler, 1571-1630) (太陽與行星運動的規律性) 與伽利略(Galileo Galilei, 1564-1642)，但最後集大成的是牛頓。大家或許都熟知運動三大定律與萬有引力定律的內容，但不一定了解它們的重要性以及對日後物理學發展的深遠影響。牛頓在物理學上的貢獻可以歸納成幾個方面：他把描述運動現象的運動學發展成動力學 (即描述運動現象與作用力的關係)，將各種運動現象的根本原因加以整合，使其適用於天體及地面上的各種運動；發明並利用微積分處理物理學的運動分析，建立動力學理論；把數學公式化的理論方法與架構運用到物理理論上，使物理理論展現出清晰、明確與嚴謹的定義，和強而有力的推演能力，並以數學關係式表達物理量間的定性與定量關係；建立動力學理論的標準模式；也是歷史上第一位提出有自然力 (萬有引力) 存在的人。所以牛頓被稱為科學史上最重要的人物之一。

牛頓不但將人類在力學上的紛雜觀念，整合為一簡潔而完整的理論，可適用於大部分的運動現象，既是古典力學的基礎，也是古典物理學的基礎，它同時也影響了古典熱力學和古典電磁學的發展。他所發展出來的理論模式是物理學的典範，牛頓本人大概也沒料到他的貢獻與影響會如此遠大。隨後發現的基本自然力，包括電磁作用力 (發生在帶電粒子與電流之間)、強作用力 (發生在原子核內及基本粒子之間)、弱作用力 (發生在核內粒子的電子放射衰變過程中)，一些相關理論也多少受其影響。我們若要了解某系統的物理現象，根據動力學的理論，其處理程序是：了解系統的組成分子與結構；掌握組成分子的質量與相互作用力的型態；正確使用運動變化法則與運動方程式；依據統計分配法則與統計資料計算結果。由運動方程式或統計分配法則與計算，我們可以獲知系統的物理現象，包括穩定狀態時的結構與性質、受外力作用時的變化情形。這種動力學模式，事實上隱藏在牛頓力學的理論體系內。構成物理學的四大力學，即是描述不同物理系統的基本動力學理論。而各種尺寸大小、作用力型態與結構不同的物理體系，也在這種基本動力學模式下，建立起本身的物理細部理論架構。 (部分文章摘自：自然科學之母—物理科學，作者：唐富欽 成功大學物理系)

第3.1節 座標系統簡介
數學上，座標系統的觀念有助於幾何問題的描述，以及對抽象概念的理解，例如一n維的向量空間，在賦予其適當之$n$個互相垂直之座標軸之後，即可視為與$\mathbb{R}^n$空間同構 (isomorphic)，$\mathbb{R}^n$是一種歐氏空間 (Euclidean Space)。數學上常用的座標系統是二維以及三維的直角 (rectangular) 座標系統，分別用來描述$\mathbb{R}^2$以及$\mathbb{R}^3$空間，是由數學家笛卡耳 (Rene Descartes, 1596-1650) 首先提出，這兩種座標系統通常也稱為卡氏座標系統 (Cartesian Coordinate Systems)。

圖3.1 平面卡氏直角座標：平面上的任一點是以序對$(a, b)$來表示

圖3.1為標準的平面卡氏直角座標，在平面上的兩條直線分別稱為$x$軸(水平線)與$y$軸(鉛直線)，其交點$O$稱為原點。以原點為準，$x$軸的右方為正，左方為負，$y$軸的上方為正，下方為負。平面上的任意一點$P$，向$x$軸與$y$軸各作垂線，假設$P$在$x$軸與$y$軸之垂足分別為$a$及$b$，則$P$點的平面座標就以序對 (ordered pair) $(a, b)$ 表示。所謂序對係指$a$, $b$是不可交換的，也就是說 $(a, b)\ne (b, a)$，而且$(a_1, b_1)= (a_2, b_2)$之充分且必要條件為$a_1=a_2$且$b_1=b_2$。如同前面所說的，$\mathbb{R}^2$是一個歐氏空間。在歐氏空間上，兩點可以決定一直線，而且兩點之間的距離可以利用畢氏定理 (Pythagoras Theorem) 求得 (如圖 3.2)，假設兩個點$A$, $B$的座標分別為$(a_1, b_1)$以及$(a_2, b_2)$，則
\[\overline{AB} = \sqrt{(a_1-a_2)^2 + (b_1 - b_2)^2}. \hspace{3cm}\mbox{(3.1-1}) \]
在日常生活中，較常使用的則是三維的直角座標，在$\mathbb{R}^3$空間中的任一點是由序對 $(a, b, c)$ 表示，同樣地$a$, $b$, $c$是不可交換的。$\mathbb{R}^2$也是歐氏空間，因此兩點亦可決定一直線，而兩點的距離可以利用畢氏定理決定：

圖3.2 平面直角座標。任一兩點可以決定一直線，其距離可以利用畢氏定理求得

圖3.3 三度空間卡氏直角座標系統：空間中任一點是由序對 $(a, b, c)$ 表示

\[\overline{AB} = \sqrt{(a_1-a_2)^2 + (b_1 - b_2)^2 + (c_1 - c_2)^2}. \hspace{3cm}\mbox{(3.1-2}) \]

另一種常用的座標軸是極座標系統，在此系統中，空間中的任意一點是用半徑以及其與直角座標軸之間的夾角來表示。如圖3.4，在平面上取一固定點為原點，向右的水平射線$\overline{OA}$，稱此射線為極軸 (polar axis)，$O$為極 (pole)，$P$為平面上任意一點，若線段$\overline{OP}$之長度為$r (r\ge 0)$，而射線$\overline{OA}$與$\overline{OP}$之夾角為$\theta$，則以$(r,\theta)$來表示$P$點之座標，此稱為$P$點之極座標。極座標的一個特點是座標的表示方式並不是唯一的，事實上，$P$點之極座標亦可以$(r, 2\pi+\theta)$表示，其中$2\pi+\theta$為$\theta$之同界角。三度空間的極座標則需要r以及兩個角度來表示，如圖3.5所示。

圖3.4 平面座標系 – 極座標

圖3.5 三度空間極座標系統

直角座標系統與極座標系統之間可透過簡單的三角函數關係，進行轉換。假設平面上一點P之直角座標與極座標分別為$(x, y)$及$(r,\theta)$，則其間關係如下 (參考圖3.4)：

\[\left\{ \begin{array}{l}
x = r\cos \theta \\
y = r\sin \theta
\end{array} \right.,\;\;\;\left\{ \begin{array}{l}
\theta = {\tan ^{ - 1}}\frac{y}{x}\\
r = \sqrt {{x^2} + {y^2}}
\end{array} \right.,\;\;\;\left\{ \begin{array}{l}
{\rm{sin}}\theta {\rm{ = }}\frac{y}{{\sqrt {{x^2} + {y^2}} }}\\
\cos \theta = \frac{x}{{\sqrt {{x^2} + {y^2}} }}
\end{array} \right..\,\,\,\,\,\, (3.1-3)\]

同理，假設三度空間上一點P之直角座標與極座標分別為$(x, y, z)$及$(r,\theta,\phi)$，則其間關係可以表示如下(參考圖3.5)：

\[\left\{ \begin{array}{l}
x = r\cos \theta \cos \phi \\
y = r\sin \theta \cos \phi \\
z = r\sin \phi
\end{array} \right.,\,\,\,\,\left\{ \begin{array}{l}
\theta = {\tan ^{ - 1}}\frac{y}{x}\\
r = \sqrt {{x^2} + {y^2} + {z^2}} \\
\phi = {\tan ^{ - 1}}\left( {\frac{z}{{\sqrt {{x^2} + {y^2}} }}} \right)
\end{array} \right.\,\,\,\,\, (3.1-4)\]

導航是什麼？GPS又是什麼碗稞，可以吃嗎？

導航 (Navigation) 簡單的說就是要知道自己所在的位置，另一個類似的名詞叫做監視 (Surveillance)，則是要知道別人所在的位置。監視所使用的工具大多為雷達 (Radar)，不在本文章所討論的範圍。

導航系統的歷史已經十分悠久，例如黃帝時代的指南車也可視為導航系統的一種；不過，事實上指南車只能指出方向，無法得知現在所在的位置。大概在十五、十六世紀左右，導航主要還是指船舶在海上航行的方向與位置導引，導航的方式係以天文導航 (celestial navigation) 為主，亦即觀測天上的恆星，以判斷載具前進的方向。這可以從導航這個字的字根 (navi) 看出來，在拉丁文中，"navis"是船舶的意思。既然是觀測天上的星星，這個時期的導航系統就是一些利於觀星的裝置，例如：Astrolabe, Sextant。

一直到20世紀初，無線電傳輸的技術成熟，導航系統才有突破性的進展；同時，由於飛機的發明，導航的方式也有重大的改變。這時期的導航設備，大多是所謂的無線電導航系統 (Radio Navigation System)。這裡我們可以先來談一談怎麼利用無線電系統來定出自己的位置，其實原理說穿了一點也不困難，簡單的說就是利用距離以及方向的資訊。無線電導航系統通常需要有參考站，設置在已知的位置上，亦即參考站的三度空間座標是已知的，該參考站持續發射一個經過特殊方式調變過的無線電信號。而在載具上 (例如飛機或者船舶) 則需要裝設無線電信號的接收裝置，這個導航設備在接收到參考站發射的無線電信號之後，由於該信號的特殊調變方式，使用者可據以計算載具與參考站之間的距離以及相對的方位。由於參考站的座標為已知，因此可以推算出載具的三度空間座標。這種定位方式稱為距離-方位定位法 (rho-theta positioning；另外幾種常見的方式稱為rho-rho positioning、theta-theta positioning以及hyperbolic positioning，將在稍後介紹)。無線電導航系統的種類非常多，例如VOR/DME、Omega、Loran-C、ILS、MLS等。我想這些名稱可能有些人一輩子也沒聽過，或許只有民航局官員、飛機駕駛員還有船長等專門的人員才會接觸到這些東西。這很正常，因為這些導航系統是給飛機以及船舶使用的，一般普羅大眾根本不需要用到這些東西。大概在10年以前，應該還沒有人會預見導航系統竟然會成為隨處可見的日常民生用品，這都要感謝GPS的發明！

在GPS (Global Positioning System，全球定位系統) 發明之前，導航一直都不是一門大眾科學，甚至在GPS發明之後的十年內，多數人還是認為：開車幹嘛用導航系統，手機裝GPS做甚麼，又不能吃。傳統的無線電導航設備價格都非常昂貴，而且體積龐大，當然不可能變成日常的消費電子產品。早期的GPS接收機體積也很大，同樣不適合於隨身攜帶。隨著半導體製程的技術越來越進步，GPS接收機的體積也越做越小，現在的技術已經可以做成單一顆IC，可以輕易的裝在手機裡面，而且價格也很便宜，現在來看，GPS接收機已經是一種消費性電子產品了。不過，我認為GPS只所以現在變得這麼普遍，有一個重要的原因應該是因為GPS所提供的服務是免費的。前面所提到的那些無線電導航設備 (VOR/DME, ILS等)，航空公司的飛機在巡航時，接收到這些參考站的信號進行定位，需要付錢給各國的民航官方單位 (因為這些參考站通常是由各個國家自行建置，由官方的機構，通常就是民航局來管理，提供導航服務，這些服務是要收費的)，收費多半非常的高昂 (航空公司應該也不在意收費多高吧，反正都會轉嫁給消費者)。GPS也可以看成是一種參考站，只不過這個參考站是架設在太空中的衛星，由衛星發射無線電信號，提供定位服務。GPS是美國國防部的財產，美國同意免費提供GPS信號給全世界使用，我認為這也是GPS迅速成為最普遍的導航設備的一個重要原因。不過，相對於傳統的無線電導航系統，GPS還是佔有許多優勢，才有辦法迅速的取代舊式的導航系統。

圖1. 利用多邊定位的方法，計算Wing Gundam的位置


GPS導航的方式係採用rho-rho positioning，也就是距離修正的定位方式 (如圖1)。若參考站A與B之位置為已知，Wing Gundam與A之距離為rho1，與B之距離為rho2，則 Wing Gundam應位於分別以A, B為圓心，rho1, rho2為半徑的兩個圓的交點上。

GPS衛星訊號提供兩種基本的觀測量 (observables)：電碼相位觀測量 (code phase observable) 以及載波相位觀測量 (carrier phase observable)。電碼追蹤 (code tracking) 結果可得到電碼相位觀測量 (code phase observable)，可用來估測衛星與接收機之間的 (瞬時) 相對距離，因為電碼相位觀測量在傳遞過程中會發生延遲，因此所得到的並非真實距離 (true range or geometric range) 而是所謂的虛擬距離 (pseudorange)。另一種觀測量則是載波相位觀測量 (carrier phase observable)，是利用接收到之衛星訊號的載波相位，以及接收機產生之參考訊號的載波相位之間的相位偏移 (carrier phase offset) 來計算衛星與接收機之間的相對距離，也可以用來計算距離的變化率，亦即都卜勒頻率 (Doppler frequency)。利用載波相位觀測量可以得到較為精確的定位結果，但是因為其訊號為週期波的關係，會有周波未定值 (cycle ambiguity) 的問題。)

訊號與系統 (I)

1.1 信號模型 (Signal Models)
常見的訊號分類方式有下列幾種：
1. 連續時間 (continuous-time) 與離散時間 (discrete-time) 信號
2. 類比 (analog) 與數位 (digital) 信號
3. 週期 (periodic) 與非週期 (aperiodic) 信號
4. 能量 (energy) 與功率 (power) 信號
5. 決定性 (deterministic )與隨機 (stochastic) 信號

簡單說明如下：
1.1.1 連續時間與離散時間信號
如同字面上的意思，若一個信號的值可以在連續的時間點上指定，稱為連續時間信號，例如：
\[x(t)=e^t,\,\,\, t\in\mathbb{R}.\]
若信號的值僅存在於離散的時間點上則稱為離散時間信號，例如：
\[u[n]=\sin n, \,\,\, n=0, 1, 2, 3,...\]
1.1.2 類比與數位信號
一般同學容易將類比信號與連續時間信號搞混，這兩種信號分類並不相同；同樣地，數位信號與離散時間信號也是不同的信號分類。若一個信號的振幅可以是任意的連續數值 (當然一般而言是在某一個範圍之內)則稱為類比信號，亦即該信號的可能數值有無限多個；反之，若該信號的振幅僅可能是有限個數值中的某一個值，則稱為數位信號。若是從函數的觀點來看，對一個信號而言，其振幅可以視為是時間的函數，亦即時間是自變數，振幅是應變數。因此，自變數是連續的為連續時間信號，自變數是離散的則為離散時間信號；應變數為連續的稱為類比信號，應變數為離散的稱為數位信號 (參考圖1，取自B.P. Lathi, "Modern Digital and Analog Communication Systems)。

圖1.1 四種不同的信號：(a)類比、連續時間信號 (b)數位、連續時間信號
(c)類比、離散時間信號 (d)數位、離散時間信號

1.1.3 週期與非週期信號
對一個信號$g(t)$而言，若存在正值$T_0$，使得$g(t)=g(t+T_0)$，對所有的$t$均成立，則稱$g(t)$為週期信號 (periodic signal)，滿足上述條件的最小正數$T_0$，稱為$g(t)$之週期。若$g(t)$不為週期信號，則稱為非週期 (aperiodic) 信號。

1.1.4 能量與功率信號
對一個信號$g(t)$而言，其能量 (Energy) 定義為：
\[E=\int^{\infty}_{-\infty}\vert g(t)\vert^2 dt\]
功率 (power) 係定義為
\[P=\lim_{T\rightarrow\infty}\frac{1}{T}\int^{T/2}_{-T/2}\vert g(t)\vert^2 dt\]
若$g(t)$的能量為有限，則稱為能量信號；反之，若$g(t)$的功率不為零，且為有限數值，則稱為功率信號。觀察上列兩個方程式可以看出來，功率係為能量的時間平均。由於時間平均的區間為無限大，因此若$g(t)$為能量信號，其功率必為零；反之，若$g(t)$的功率為有限，則其能量必為無限大。換句話說，一個信號不可能同時是能量信號又是功率信號；但有可能既不是能量信號也不是功率信號，例如斜坡信號 (ramp signal)。

1.1.5 決定性與隨機信號
若一個信號可以利用物理或者數學方式完整描述，或者可以利用時間的函數來表示，在任何時間點其函數值是可確定的，稱為決定性信號。反之，若在任一時刻，其數值僅能用機率的方式來描述，例如期望值、標準差等，則稱為隨機信號。在通訊原理中，一般說來系統中的雜訊都是以隨機信號的方式來表示。

Example: Find $P$ and $E$ for the signal $x(t)=e^{-2t}u(t)$, where $u(t)$ is the unit step function.
Solution:
\[ E=\int^{\infty}_{-\infty}\vert e^{-2t}u(t)\vert^2 dt = \int^{\infty}_0 e^{-4t}dt = \frac{1}{4} \]
\[ P=\lim_{T\rightarrow\infty}\frac{1}{2T}\int^T_{-T}\vert e^{-2t}u(t)\vert^2 dt = \lim_{T\rightarrow\infty}\frac{1}{2T}\int^T_0 e^{-4t} dt = \lim_{T\rightarrow\infty}\frac{1-e^{-4T}}{8T} = 0\]
Therefore, $x(t)$ is an energy signal.

Stationary Random Process

在工程上，特別是通訊工程，我們經常碰到隨機程序；而為了處理上的方便，通常會假設隨機程序具有stationary，亦即穩態的特性。何謂穩態，簡單的說就是我們所作的實驗不會因為作實驗的時間不同，而得到不同的結論。例如擲骰子，早上擲和晚上擲，擲出任何一面的機率都是六分之一。當然這樣的假設不一定永遠成立，例如GPS定位，一般而言晚上作實驗定位誤差會比早上作實驗來的小。

以數學的方式來描述穩態特性，假設$X(t+1), X(t+2), ..., X(t+n)$為隨機程序$X$在$n$個時間點上的採樣，其聯合機率密度函數為$f(X(t+1), X(t+2),..., X(t+n))$，則$X(t)$為stationary的條件為：

\[f(X(t+1), X(t+2),..., X(t+n)) \\ \,\,\,\,\,\,\,\, = f(X(t+1+k), X(t+2+k),..., X(t+n+k)) \\ \,\,\,\, \mbox{對所有的$k$均成立 (Equation 1)} \]
(PS: 這樣的寫法並不十分正確，此處僅希望說明stationary的概念，我們就先不管數學嚴不嚴僅吧！)

也就是說，隨機程序的聯合機率密度函數不會隨著時間改變。在線性系統理論中，這有點類似非時變 (time-invariant) 的概念，不過對於一個隨機的東西，稱它是非時變好像有點奇怪，因此比較好的名稱可以說它是shift-invariant。不過，stationary random process跟線性非時變 (Linear Time-Invariant, LTI) 系統有密切的關係：一個高斯白雜訊通過一個線性非時變系統，其輸出必為一個stationary random process。

Equation 1這樣的條件，看起來似乎沒什麼問題，但實際上確不怎麼好用。機率密度函數可以完整描述一個隨機變數，但是實際上並不容易求得隨機變數的機率密度函數，更重要的是在工程應用上，大部份的情況下我們並不需要知道隨機變數或者隨機程序完整的機率密度函數，通常我們僅需要知道它的一階以及二階統計特性 (first- and second-order moment)，亦即所謂的期望值與變異量 (或者叫做協方差)。這就好像電子學裡面的小訊號模型，當系統的交流信號值很小時，我們僅需考慮系統的一階微分項，此時系統的輸入與輸出的關係僅是一個單純的線性映射，這個線性映射就是Two-Port Network，也就是小信號模型。再回到隨機程序的問題，shift-invariant的性質，在一階以及二階統計特性所體現出來的事實就是，期望值必須為常數，以及自相關函數 (autocorrelation) 必須僅與時間差有關 (shift-invariant在三階以上的moment所體現出來的特性則較少在文獻上討論)。這就衍伸出wide-sense stationary (WSS) 的概念，我們稱滿足 (1) 期望值為常數 (2) 自相關函數僅與時間差有關 這兩個特性的隨機程序為 WSS隨機程序。有時為了區分，我們會稱前面介紹的stationary性質為strict-sense stationary (SSS)。很明顯的SSS隨機程序必為WSS隨機程序，反之則未必成立。

這邊需要澄清一下，WSS僅需滿足上面(1)(2)的性質，因此其三階以上統計特性shift-invariant的性質未必滿足。同時，(1)(2)兩個性質是shift-invariant所造成的結果，(1)和(2)之間未必存在因果關係，換句話說，期望值為常數未必自相關函數就會僅與時間差有關，反之，自相關函數僅與時間差有關，未必期望值就會是常數。

接下來我想談一下Wiener Filter以及Kalman Filter的一些觀念。(待續)

LaTeX works on blogger!

\[
\begin{align}
\nabla \times \vec{\mathbf{B}} -\, \frac1c\, \frac{\partial\vec{\mathbf{E}}}
{\partial t} & = \frac{4\pi}{c}\vec{\mathbf{j}} \\
\nabla \cdot \vec{\mathbf{E}} & = 4 \pi \rho \\
\nabla \times \vec{\mathbf{E}}\, +\, \frac1c\, \frac{\partial\vec{\mathbf{B}}}{\partial t} & = \vec{\mathbf{0}} \\
\nabla \cdot \vec{\mathbf{B}} & = 0
\end{align}
\]