行政院國家科學委員會專題研究計畫成果期中報告
一個表達三維城市表面模型的小波最小二乘法(1/2)
A Wavelets-Based Least-Squares Algorithm for 3D-City Modeling (1/2)
計畫編號：NSC 89-2211-E-006-082
執行期限：88 年 8 月 1 日至 89 年 7 月 31 日
主 持 人：蔡展榮 執行機構及單位名稱：國立成功大學測量工程學系
一、中文摘要
在第一年度裡，本研究已對若干基本問題做 一詳細研究，並獲致若干具體結論，陸續在國內 外的學術研討會與論文期刊發表研究成果。例 如，運用各式模擬以及實測的資料，來對以小波 表達一個具有下列三項特性的城市幾何表面之 精度、收斂性與可行性進行分析：À訊號各處的 連續度可不相同、Á包含（似）斷點或（似）斷 線之局部訊號、Â多層解析度之訊號。研究結果 顯示，小波法確實能對含似斷線和似斷點之幾何 表面作整體精確的表達，而不需分區地描述。小 波表達之近似函數的精度也確實會隨著解析度 的變細而提高，其於局部似斷點和似斷線的表達 上，吉布斯現象的擺振範圍也會隨著解析度的變 細（即：解析度指標j值變大）而以2
−j倍迅速 縮小。而且近似函數於似斷點處的局部最大誤差 量與局部高差之比值是呈一常數。具碎形特性之 Daubechies 三階小波對平滑、崎嶇不規則與含似 斷點之訊號都能精確地表達，且相較於 Haar 小 波與 Meyer 小波則有較大之適用性。小波對具有 多層解析度性質的真實表面之表達來得比傳統 傅立葉表達法更精確且更有效率。關鍵詞
：三維城市表面模型、小波、最小二乘法Abstr act
In the first year, some fundamental points are studied. The results are consecutively published in some domestic and foreign conferences and journals. For instance, some simulated and real test data are used to do a van study. To exam the characteristics of wavelet-based signal representa-tion, one focuses mainly on those
signals having the following properties: • various degree of continuity on different signal points, ‚ having local signals with (pseudo) break points or –lines, ƒ multi-resolution signals. It is concluded that a wavelet based approximation can describe an entire geometrical signal with local (pseudo) break points and/or –lines, where piecewise representation is not needed. The accuracy of such an approximation becomes better if a finer resolution (i.e. a larger j-value) is adopted. The influence radius of the Gibbs phenomenon becomes smaller in a rate of
2
−j as the j-value becomes larger. The ratio of the maximal approximation error near a pseudo break point to the signal height difference at that break point is a constant. On the other hand, the Daubechies wavelet of 3rd order is the best basis comparing to the other representative bases, namely the Haar wavelet, the Meyer wavelet and the Fourier basis. It is very suitable to describe different kinds of signals, e.g. a smooth signal, a rugged one and a signal with local break points. Generally speaking, the wavelet representation is more accurate and efficient than the conventional Fourier method.Keywor ds
: 3D city surface modeling, wavelets, Least Squares Algorithm二、緣由與目的
本文的研究目的是在•研究一個有效表達 三維數值城市表面模型的小波最小二乘法，以其 能夠提出一套具有實用價值的數學模式與計算 方法、‚建立以三維城市表面上離散分佈的點群 來表達城市模型的數學模式，並以最小二乘原理 來求定相關的模式參數、ƒ做為後續的自動化重 建三維城市模型研究之先驅研究。三、數學模式
在本文中，將分別以符號 Z、R 表示全部整 數、實數之集合。真實世界的許多訊號f是典型 的 「 多 層 解 析 度 訊 號 (multi-resolution signals)」，訊號各處具有粗細不同的「解析 度」，例如，一般的城市幾何表面就是一例，有 些地方是大面積之高樓大廈區，有些地方是小面 積之低矮房子區，更有些地方是平坦之廣場… 等 等，它們常是一種碎形、各處連續度可互不相同 的多層解析度訊號。 當一個實數訊號函數 f 滿足式(1)之條件 時，則 f 可分解出不同解析度下的近似訊號
f
A
j ，而A
jf
可依式(2)所示的一個線性式來表 示之[4,5]：∫
+∞( )
∞ − f x dx 2 < ∞ (1) Af(x) a jk(x) z k jk j∑
φ ∈ = (2) 其中，j 是一個表達解析度的指標，∀j∈Z，j 值愈大，表示解析度愈細；φjk(x)是一個已知的 基 底 函 數 ， 稱 為 「 尺 度 函 數 (scaling function)」，它是由「小波父函數φ(x)(father wavelet)」依式(3)作不同的平移和尺度調整而 成，ajk = f,φjk =∫
f(x)φjk(x)dx +∞ ∞ − ，稱之為「尺 度係數(scaling coefficients)」[1]。 (x) 2j (2jx k) jk = φ − φ ，∀j,k∈Z (3) 另外，f 也可分解出不同解析度下的細節訊號f
D
j ，如式(4)所示的一個線性式來表示之： D f(x) w jk(x) z k jk j∑
ψ ∈ = (4) 其中，ψjk(x)是另一種已知的基底函數，稱為「小 波函數」，它是由「小波母函數ψ(x)(mother wavelet)」依式(5)作不同的平移和尺度調整而 成，wjk = f,ψjk =∫
f(x)ψjk(x)dx +∞ ∞ − ，稱之為「小 波係數(wavelet coefficients)」[5]。 ) 2 ( 2 ) (x j jx k jk = ψ − ψ ，∀j,k∈Z (5) 根 據 多 層 解 析 度 分 析 (multiresolution analysis, MRA)之理論，前述的訊號函數f可表 示成： f(x) D f w jk(x) j k jk j j∑∑
ψ∑
=+∞−∞ = +∞=−∞ +∞=−∞ = (6) 而 上 式 就 是 所 謂 的 為 「 小 波 級 數 (wavelet series)」[4]。另一方面，根據 MRA 理論，式(6) 也等於：∑
+∞ = + = 0 0 ( ) ) ( j j j j f x D f A x f (7) 上列各式就是採用小波來表達（城市）幾何表面 剖線函數之基礎理論公式。四、小波表達幾何表面函數之方法（一）
選定初始解析度 j 求定平移參數k之值域 [kmin,kmax] 以數值積分法估算全 部尺度係數ajk 計算並輸出近似函數 Ajf max∣Aj+1 f (x)-Aj f (x)∣ <ε? 是 結束 輸入全部參考 點之參考值 選擇小波之種類 否 j=j+1 給定收斂臨界值ε 圖 1. 小波表達幾何表面函數之方法（一）的演 算流程 在本文中，採用數值積分法來求定表達式(2) 裡的各個尺度係數ajk，其中，平移參數 k 的值 域可根據欲分析的訊號 f 的分佈區間[x
1,x
n] 依照式(8)來求得： 2 x1 R1 k 2 xn R2 j j − ≤ ≤ − (8) 其中，[R1,R2]表示小波父函數φ(x)的承載區 (support)。此一小波表達幾何表面函數之方法 （一）的演算流程如圖 1 所示。五、精度影響因子之分析
以城市幾何表面的一維剖線訊號為基礎，針 對訊號本身之斷點、坡度、與小波近似法中的解 析度指標 j、基底函數之種類、以及訊號的類型 等五種主要影響剖線表達精度之因子作分析，並 由分析中探討其對所表達的城市幾何表面剖線 精度的影響與小波表達城市幾何表面函數的方 法之收斂性與經濟效益。 5－1 斷點與坡度 在本文的研究中，調整階梯訊號裡的斷點處 上下兩端點之水平間距△x，分析其對近似函數 的影響，得到下列幾點具體的結論： •（似）斷點處上下兩端點之吉布斯現象(Gibbs phenomenon)所產 生的近似訊號A
jf
擺振範圍 為： [x1-L j 2 , x2+L j 2 ] (9) 其中，符號 x1、x2 是（似）斷點處上下兩端點 之平面坐標，L為小波父函數φ(x)主要非零值分 佈區之寬度。 ‚（似）斷點處之A
jf
擺振範圍會隨著j的提高， 而以 −j 2 倍迅速縮小。 ƒ當△x＞φjk(x) 之主非零值分佈區寬度，則f
A
j 擺振現象產生之最大誤差量相較於高差便 顯得不顯著。 5－2 解析度 由模擬一個已知的城市剖線訊號（圖 2）來 分析解析度（指標 j）對近似函數Ajf的精度與 收斂性之影響。由實驗結果得到近似函數Ajf 與 原函數f之間的真誤差的均方根值與最大誤差量 和解析度指標j的關係曲線，如圖 3、4 所示。 由本節相關的實驗結果分析得出下列幾項 結論： • 近似函數Ajf 確實會隨著解析度變細而迅速 收斂。 ‚ 較低頻之局部訊號以較低解析度之近似函數 就能做很好的表達；而較高頻之局部訊號則需較 高解析度之近似函數才能作較精確的表達。 ƒ 較大之訊號表達誤差均發生於似斷點處。 „ 近似函數於似斷點處的局部最大誤差量與該 處的訊號高差之比值呈一常數，以j=4 之近似函 數為例，局部最大誤差量與高差之比值為 0.03。 … 尺度係數ajk之個數會隨著解析度指標j值的 增大而以 2 的等比級數遞增。 0 5 10 15 20 25 30 35 40 -6 -4 -2 0 2 4 6 8 10 12 x f( x ) 圖 2. 一個已知城市剖線訊號f 0 1 2 3 4 5 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 j 圖 3. Ajf 與真值f之均方根差 0 1 2 3 4 5 0 1 2 3 4 5 j 圖 4. Ajf 之最大誤差量 5－3 基底函數 為了探討不同小波基底之近似函數表達不 同性質訊號的精度影響，吾人分別以圖 5 所示的 三種具代表性的 Haar 小波（含零階不連續）、 Daubechies 三階小波（連續度α為 1.0878≦α ＜2）與 Meyer 小波（各處皆為無窮階連續）為基底來表達含似斷點之城市剖線（圖 2）、崎嶇 不規則與平滑之 Sin 函數訊號（圖 6）。圖 7 顯 示三種小波表達三種訊號的均方根差與最大誤 差量之關係曲線。 -0.5 0 0.5 1 1.5 -1.5 -1 -0.5 0 0.5 1 1.5 x a m p lit u d e Haar 尺度函數 -1 0 1 2 3 4 5 6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 x a m p lit u d e Daubechies 不 對 稱 尺 度 函 數 N=3 -8 -6 -4 -2 0 2 4 6 8 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1 1.2 x a m p litu d e Meyer 尺度 函數 圖 5. Haar
φ
(x
)（上）、Daubechies 三階φ
(x
) （中）與 Meyerφ
(x
)（下） 0 2 4 6 8 10 12 14 16 18 20 0 2 4 6 8 10 12 14 16 0 2 4 6 8 10 12 14 16 18 20 -8 -6 -4 -2 0 2 4 6 8 圖 6. 崎嶇不規則訊號（上）與 Sin 訊號函數（下） 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 j RM S E Haar Daub.3 M ey er 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 0 1 2 3 4 5 6 7 8 9 j M a x E rr o r Haar Daub.3 M ey er （a）三種小波對城市剖線表達之均方根差(上) 與最大誤差量(下)之比較 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 j R M S E Haar Daub. 3 M ey er 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 0 1 2 3 4 5 6 7 8 9 j M a x E rr o r Haar Daub. 3 M ey er （b）三種小波對崎嶇不規則訊號表達之均方根 差(上)與最大誤差量(下)之比較 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 j R M S E Haar Daub. 3 M ey er0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 0 1 2 3 4 5 6 7 8 9 j M a x E rr o r Haar Daub. 3 M ey er （c）三種小波對Sin函數訊號表達之均方根差 (上)與最大誤差量(下)之比較 圖 7. 三種小波表達三種訊號的均方根差與最大 誤差量之比較 由研究結果得知： • Haar小波對一般（不含斷點）的連續訊號皆 會產生鋸齒狀之表達函數，且其精度比其它兩者 來的差。 ‚具碎形特性之 Daubechies 三階小波對三種訊 號的表達精度均比其餘兩者來的佳，其適用性也 來的比其他兩者廣。
六、小波表達幾何表面函數之方法（二）
6-1 表達式 為改善方法（一）的計算效益而改用引自式 (7)結合小波與尺度函數之表達式，如下：∑
− = + = 1 0 0 ( ) ( ) ) ( J j j j j J f x A f x D f x A (10) 其中，A0f(x) a 0 j0k(x) z k k j j∑
φ ∈ = 且 ) ( ) (x w x f D jk z k jk j∑
ψ ∈ = 。式(10)之特點如下： •能以較低解析度之近似函數Aj 0f 來表達訊號 中的低頻分量，而以高解析度之細節函數Djf 來 表達訊號中的高頻分量，使得表達可以更有效 率。 ‚表達式中的小波係數wjk可將相鄰解析度間的 尺度係數ajk、aj+1,k之相關性加以抽離，且小波 係數只於訊號中的顯著高頻處有值，其餘於低頻 處則近乎等於零，因此，只需記錄顯著非零值即 可，這使得記錄量大幅降低。 ƒ當近似函數Ajf 逐層提高解析度時，不需如方 法（一）的尺度函數表達式，每個解析度的近似 函數需重新計算一次，而只需計算所需解析度的 顯著細節函數即可。 小波表達幾何表面函數之方法（二）的演算 流程如圖 8 所示。 輸入全部參考 點之參考值 選定使用之正交小波 基底函數 選定初始解析度j0 以數值積分法求得 aj0,k，wj0,k ，V k 計算於j之臨界值ξ ，並 過濾出顯著之小波係數 值wj,k，且記錄其位置 由顯著小波係數值 wj,k之位置來計算wj+1,k 是 結束 max∣Aj+1 f (x)-Aj f (x)∣ < ε? 否 j=j+1 給定收斂臨界值ε j=j0 計算並輸出近似函數 Aj+1 f (x) 計算並輸出近似函數 Aj f (x) 圖 8. 小波表達幾何表面函數之方法（二）的演 算流程 6-2 表達式修改前後之精度與效益 以一維城市剖線（圖 2）之近似函數A4f為 例，表達式修改前（方法一）與表達式修改後（方 法二）之比較顯示，修改後之表達式幾乎不影響 近似函數之精度，而係數個數（累加項數）卻減 少了約 81﹪。七、實際城市表面函數之測試成果與分析
以二維之實際城市表面資料來測試以小波 表達二維城市幾何表面函數之精度與收斂性。而本研究所採用的實驗區是位於新莊地區，圖 9 與 圖 10 分別是實驗區之 CAD 圖與攝影影像。採用 解析製圖儀（Leica BC3）量測的建物屋角點、 屋脊點、DTM 與位相關係線性內插得實驗區的表 面模型，如圖 11 所示，並把此表面視為一個已 知的城市幾何表面函數來做測試。 圖 12 是以方法（二）所得到的實驗區近似 表面函數Ajf（基底為 Daubechies 三階φjk(x)）。 由圖 12 可知j從 1 至 4 之近似表面函數，其幾 何形狀是漸漸地逼近原始的已知表面（圖 11 左）。 1B 3R 4R 3R 4R 3R 4R 1B 圖 9. 實驗區之 CAD 圖 圖 10. 實驗區之攝影影像 圖 11. 實驗區幾何表面模型 f A1 A2f f A3
A4f 圖 12. 近似表面函數Ajf 圖 13 是以正射影像來與j 從 0 至 4 的近似 表面函數誤差量的灰階影像做一對照比較。圖中 色調愈黑表示誤差愈大，由圖 13 可以發現誤差 較大之地方均發生於斷線處，且斷線處的高差愈 大相對地誤差也就愈大，而平滑之區域，如左下 之農田區，於j=0 時其誤差就已不明顯。由圖 13 與圖 14 也可知，隨著解析度變細，誤差量也確 實漸漸地變小。 (a) (b) (c) (d) (e) (f) 圖 13. 正射影像（a）與A0f∼A4f 的誤差量 灰階影像（b）∼（f）之比對 0 0.5 1 1.5 2 2.5 3 3.5 4 0 0.2 0.4 0.6 0.8 1 1.2 1.4 j
00
.
5
3
4
5
5
11
.
5
22
.
5
33
.
4
6
7
8
9
1
01
11
2
j
M
a
x
e
r
r
o
r
圖 14. 均方差曲線（左）與最大誤差量曲線（右）八、結論 由研究成果可歸納出小波表達城市幾何表 面函數的幾項具體結論： •當似斷點處左右兩端的兩參考點之平面位置 間距△x 大於尺度函數 之主要非零值範圍 之寬度時，近似函數於似斷點處的吉布斯現象所 產生的最大誤差量相較於高差便會變得相當不 顯著。 ‚似斷點處上下兩端點附近之近似函數擺振範 圍為[x1 - x2+, ]，其中，x1 、x2 為似斷 點處左右兩端點之平面位置坐標，L 為小波父函 數主要非零值分佈的定義域區間之寬度， 為 尺度函數 主要非零值的定義域區間之寬 度。由此一結論可知，似斷點處的吉布斯擺振範 圍會隨著解析度變細（即：解析度指標 j 值變 大），而以 倍的變化速度迅速縮小。 ƒ由小波表達虛擬城市剖線與真實城市表面的 實驗結果可知，小波表達法確實會隨著解析度變 細（即 j 值增大）而提高其表達函數之精度。也 從而證實小波確實能將各處連續度（可）不一的 訊號作一整體的表達，而不需分區分段地描述。 „於第六章二維之實驗中，於實驗區中的平緩地 區，如：左下之農田區與建物之屋頂面，於低解 析度（j=0 ）時，就可得到很好的表達，而似斷 線處之訊號是隨著 j 值從 0 增到 4 才漸漸地收 斂。因此，當解析度愈低（j 值愈小），小波近 似函數所表達的近似訊號是一個愈低頻的局部 趨勢訊號，而當 j 值提高時，小波近似函數所能 表達的訊號低頻率域也就愈寬，愈能表達起伏急 遽之訊號。而這也造成了於訊號中的每個區域都 具有不同大小的誤差量，於愈平滑（低頻）的區 域，則其誤差量愈小，反之則愈大。 …當以同一解析度下之近似函數表達城市剖線 時，近似函數在似斷點處的局部最大誤差量與高 差之比值是幾乎呈一常數。 †Haar 小波能對一樣含有零階不連續之階梯狀 函數作完全精確之表達，但對其它各式訊號皆會 產生鋸齒狀之表達函數，換言之，Haar 小波不適 合對平滑訊號作高精度之表達。而 Meyer 小波對 於平滑、含似斷點兩種訊號的表達精度，雖與 D a u b e c h i e s 三階小波差異不大，但對於崎嶇不規 則的破碎訊號而言，就來的比 Daubechies 三階 小波略差。而 Daubechies 三階小波，雖是一個 碎形小波，但它對於平滑、崎嶇不規則、含似斷 點三種具代表性的訊號都能作精確之表達，另 外，具碎形性質的 Daubechies 三階小波一樣可 以 精 確 表 達 一 個 非 常 平 滑 的 函 數 ， 因 此 ， D a u b e c h i e s 三階小波的適用性來的比其他兩種 小波來的廣且佳。 ‡小波表達法確實比傳統的傅立葉表達法更有 效率而精確。由於小波基底於定義區間外之函數 值迅速收斂為零或近乎等於零之特性，而能根據 訊號中之局部特性以適宜解析度之小波基底來 表達訊號，對於訊號中的局部低頻分量，是以較 粗解析度之近似函數來表達，訊號中的局部高頻 分量，則用高解析度之細節函數來表達，而這正 是小波的多層解析度表達與時頻分析能力優於 傅立葉之所在。傅立葉基底是無限承載的單頻週 期函數 Sin(x) 、Cos(x) ，它無法針對訊號中的局 部特性以不同之合宜解析度來表達，而只能對訊 號各處均以相同的解析度來表達訊號，無法表達 各處連續度與時頻性質（可）互不相同、甚至同 時具有碎形維度的訊號，且與具有多頻、碎形特 性的小波基底相較，Sin(x) 、Cos(x) 單頻、均勻 平滑的特性，使得其對真實地表與物表幾何面函 數的表達顯得較不適用。 目前現有的幾何表面表達法皆以離散方式 來進行，本研究嘗試分析以小波建立一個連續函 數來描述含（似）斷點和（似）斷線的訊號之可 行性，經由本研究於真值已知之下所得的上述各 項結論，可以得知小波表達法的多層解析度表 達、極佳之時頻定位分析能力以及碎形描述的能 力，對於一樣具有多層解析度性質之各種「真實」 地 表 面 和 物 表 面 可 作 有效率而精確之整體表 達，證實以連續函數也可達到良好的表達，且對 局部含有（似）斷點和（似）斷線與局部平滑的 訊號一樣可同時作良好的描述。而以本文所列舉
的城市幾何表面為例，顯示此法可行。 本文所提出之演算法是根據於真實幾何表 面已知之下，以數值積分法來求得小波級數所需 之係數，然而在實務中，真實表面往往是未知 的，因此，未來如何於真實表面未知之下，對其 作最佳之取樣，然後根據取樣資料以小波表達式 建立最小二乘原理中之函數模式，配合「由粗而 細」的策略，以求解出所需之小波係數而重建最 佳表面函數，是未來後續之一項研究重點。另 外，在這一二年裡，國外陸續有了新的研究成 果，例如，[6,7,8,9,10] 。為求迎頭趕上，因此， 如何建立更佳的數學模式與演算法來處理一個 三維的多值函數(multi -valued i)noutcnf 、並 真正表達一個三維的城市表面，也是未來的一項 研究重點。
參考文獻
1、單維彰, 1998 ：“凌波初步”,全華圖書公司. 2、曾義星, 1998 ：“地理資訊系統資料處理”, 上課講義,國立成功大學.3、Kaiser, G., 1994：“A Friendly Guide to Wavelets”, Birkhaeuser.
4 、 Chui, C.K., 1992a ： “ An Introduction to Wavelets：Wavelet Analysis and its Application”, Vol. 1, Academic Press, Inc.
5 、 Daubechies, I., 1994 ： “ Ten Lectures on Wavelets”, Third Printing, Society for Industrial and Applied Mathematics, Philadelphia, Pennsylvania..
6. Foestner, W., and Guelch, E., 1999, “Automatic orientation and recognition in highly structured scenes”, ISPRS Journal of Photogrammetry and Remote Sensing 54, pp. 23-34.
7. Davis, T.J., 1999, “Fast Decomposition of Digital Curves into Polygons Using the Haar Transform”, IEEE Transactions on Pattern Analysis and Machine Intelligence, Vol. 21, No. 8, pp. 786-790.
8. Schroeder, P. and Sweldens, W., 1999,
“Spherical Wavelets: Efficiently Representing Functions on the Sphere”, University of South Carolina.
9. Pastor, L. , and Rodriguez, A., 1999, “Surface Approximation of 3D Objects from Irregularly Sampled Clouds of 3D Points Using Spherical Wavelets”, IEEE, pp. 70-75.
10. Schlueter, M., 1999, “Von der 2 1/2D- zur 3D-Flaechenmodellierung fuer die photogrammetrische Rekonstruktion im Objektraum”, DGK, Reihe C, Heft Nr. 506, Munich.
